NASA TECHNICAL 
MEMORANDUM 


NASA TM X -53389 

OCTOBER 15, 1965 


©N 

00 

CO 

CO 

1-0 

» 

X 


C/3 

z 




GPO PRICE $ 

CFSTI PRICE(S) $ 

Hard copy (HC) ^ 30 

Microfiche (MF) /. 2) Q 

ft 653 July 65 


AERO-ASTRODYNAMICS 
RESEARCH REVIEW NO. 3 




AERO-ASTRODYNAMICS LABORATORY 
RESEARCH AND DEVELOPMENT OPERATIONS 
GEORGE C. MARSHALL SPACE FLIGHT CENTER 

HUNTSVILLE, ALABAMA 



NASA - GEORGE C. MARSHALL SPACE FLIGHT CENTER 


TECHNICAL MEMORANDUM X- 53389 



RESEARCH REVIEW NUMBER THREE 
January 1, 1965 -June 30, 1965 


RESEARCH AND DEVELOPMENT OPERATIONS 
AERO-ASTRODYNAMICS LABORATORY 



ACKNOWLEDGEMENTS 


The articles for this review were contributed by var- 
ious engineers and physicists of the Aero-Astrody- 
namics Laboratory, reviewed and compiled by William 
D. Murphree, and edited by Sarah Hightower. 

Grateful acknowledgement is given to the Special 
Publications Section, Scientific and Technical Infor- 
mation Branch, Management Services Office, MSFC ; 
for preparing the review. 


ii 



PREFACE 


The topics discussed in this third Aero-Astrody- 
uamics Research Review cover a variety of subjects. 
Included are Acoustics, Aerodynamics, Base Heating, 
Control, Guidance, Instrumentation, Structural Dy- 
namics, and Orbit Theory. Although many subjects 
are discussed, all areas of research which are con- 
ducted by this laboratory are not treated. Other 
subjects will be discussed in forthcoming reviews. 

It is hoped that these reviews will be interesting 
and helpful to other organizations engaged in space 
flight research and related efforts. Criticisms of this 
review and discussions concerning individual papers 
with respective authors are invited, 

E. D. GEISSLER 

Director, Aero-Astrodynamics Laboratory 
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I. ACOUSTICS 



SOUND PRESSURE ESTIMATIONS FROM A POINT DIRECTIONAL ACOUSTIC SOURCE RADIATING 

IN AN INHOMOGENEOUS MEDIUM 



by 

Jess H, Jones 


Symbol 


Definition 


This report concerns a general expression that 
can be used to make estimates of sound pressures 
received from a rocket «xhaust that radiates in an in- 
homogeneous medium. The concepts concerning 
acoustic wave propagation are combined with those 
having to do with a broad frequency band, directional 
acoustic source, etc.; in brief, rocket exhaust flow. 
The ray acoustics approach is used in this paper; there 
is no attempt to elaborate on the failure of the race 
trace method to obtain pressure levels for special con- 
ditions. The ray trace technique is used only to de- 
scribe the relationship between horizontal source-to- 
re ceiver distance and the angle of emergence of the 
sound ray. 
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Symbol 


LIST OF SYMBOLS 


Definition 


P 2 /cps 


ref 


ASL 


Space-time mean squared sound pres- 
sure per unit bandwidth 

' oot mean squared reference sound 
pressure = 0,00002 newton/ ( meter) 2 

Power spectrum level 

r Source -to -receiver distance 

rv * 

S = Pvcps Time mean squared sound pressure per 
unit bandwiuth 

SPL Space average sound pressuie level 

SPL Sound pressure le if el 

r Radius of source 

o 

X Horizontal source -to-receiver distance 


A(f) 

B 

D.l. 

x? 

"A 

f 

Af 

H 


Power radiated per unit bandwidth 
from an equivalent monopole source 

Area 

Directivity index 
Attenuation factor 
Frequency 
Bandwidth 

Height of ray tube cross section above 
source elevation 


a 

P 

y 

6 

£ 

0 

1 


Azimuth angle of sound ray 

Directivity angle of sound source 

v 

Attenuation function 
Azimuth angle of exhaur* axis 
Elevation angle of exhaust axis 
Elevation angle of sound ray 
Attenuation coefficient in nepers/meter 
Attenuation coefficient in decibels/meter 


I 

Intensity 

m f 

Meteorological factor 

OB SPL 

Octave band sound pressure level 

OB AL 

Octave band power level 

PSL 

Pressure spectrum level 


pc Specific acc* ^iic impedance of the 

atmosphere 

0 Direc t; vity function 

0 Meteorologic.il function 

^ Summation of all rays arriving at a 

a mmon location 
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SUBSCRIPTS 


SECTION II. HOMOGENEOUS SOUND PRESSURE 
EQUATION 


* Homogeneous conditions 

0 Conditions at source 

li Conditions at source for a homogeneous 

medium 

R Conditions at receiver 

SUB-SUBSCRIPT 

1 Rays leaving the source arriving at a 

common location 


SECTION I. INTRODUCTION 


The purpose of this report is to obtain a general 
expression which can be used to make engineering es- 
timates of the far-field sound pressures received from 
a rocket exhaust acoustic source radiating in an in- 
homogeneous medium [ 1,2] . This will be accomplished 
by combining the concepts concerning acoustic wave 
propagation characteristics [3, 4, 5,6,7] with those 
concerning a broad frequency band, point, directional 
acoustic source, i.e. , rocket exhaust flow [8,9]. 

This discussion will utilize the method of ray 
acoustics and will not attempt to elaborate on the fail- 
ure of the race trace method in arriving at meaningful 
pressure levels for special conditions; i.e. , unlimited 
value of the sound pressure at a focus and scattering 
of sound into a shadow zone, etc. These areas are 
important but they are considered as being beyond the 
scope of the original intent of this report. The ray 
trace technique will be used only to describe the re- 
lationship between the horizontal source -to-receiver 
distance, X c , and the angle of emergence of the sound 

o 

ray, 0 , i.e. , X =f(0 ). 

O O o 


Under the ideal conditions of a homogeneous at- 
mosphere, the expression relating the sound pressure 
to the acoustic source can be derived with little diffi- 
culty [9,10]. The resultant equation presented here 
for a point, directional source radiating in a half space 
(including atmospheric absorption) if 


pc A(f) <p(fJ3) 

9 ] 

S*(f,j3,r) = 

27 r r 2 y(f,r) 


( 1 ) 


The directivity function <£(f,/3) describes the 
angular deviation of the radiation pattern of the actual 
rocket source from an equivalent monopole acoustic 
source A(f) . This function modifies the simple mono- 
pole source radiation pattern such that the result de- 
scribes the complex radiation pattern of the rocket 
exhaust source. The directivity function is defined as 

P 2 (f,/3) r 

4>(f,(3) = — ~ 

P 2 (f) 

r l 


IV. 

where P z (f,/3)- is assumed to be symmetrical about 
r l 

the rocket exhaust stream and for a given source fre- 
quency is a function only of the angle 0(Fig, 1 ). 
This is evaluated at a far-field distance r 1 which is 
normally kept constant with respect to 0. P z (f) r is 

the space-time mean squared pressure averaged over 
all space on the surface of a sphere of radius r*and 
consequently is independent of angular orientation. 
The space average mean squared pressure can be 
thought of as that due to an equivalent monopole source. 
For instance, for a monopole source the radiated en- 
ergy, F 2 (f,/3) , is independent of p and the directivity 
function is therefore unity. 


For the acoustic source in question, i.e., the 
rocket exhaust flow, the source-to-receiver distance 
is assumed to be large such that the locally distributed 
rocket exhausl (acoustic) source appears as a point. 
Further, the source will be considered as fixed in 
space and will exhibit highly directional pressure dis- 
tribution characteristics. 


The attenuatioi function y(f,r) is normally ex- 
pressed in exponential form and is defined here as 

y (f, r) *e 2 ^^ r , where £(f) is the attenuation coeffi- 
expressed in nepers/meters, which is assumed to be 
a function of frequency only and r is the source-to- 
receiver distance. 
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Figure 1 delineates the coordinate system used 
and also provides an expression relating the direc- 
tivity angle to the indicated angular functions. 



conditions are achieved when ip = 1. 

To evaluate ip § the sound intensity under inhomo- 
geneous conditions needs to be investigated. This is 
accomplished by utilizing the concept of ray acoustics. 
Application of this ray trace procedure is valid if it is 
assumed that the changes in the properties of the at- 
mosphere are small within a wavelength of the propa- 
gating acoustic energy. If this assumption is not sat- 
isfied, which is the case at a focus, a modification to 
the ray trace method is required. 

By studying the deformation of a ray tube (Fig. 2) 


cos IS ~ cosf cos# cos{a-6) - sin0 sin4 
o o 

FIGURE 1. COORDINATE SYSTEM 

Equation 1 provides a means of obtaining estimates 
of the mean squared sound pressure under homogen- 
eous conditions. However, it is known that when the 
atmosphere departs from this ideal condition the mean 
squared sound pressures received at a remote location 
can vary greatly from the estimates computed by using 
equation 1. Consequently, a more general form of this 
expression is highly desirable. 


y 



SECTION in. EVALUATION OF THE METERO- 
LOGICAL FUNCTION 


FIGURE 2. RAY TUBE 


It would be highly advantageous to compute inde- 
pendently the sound pressure or sound intensity due to 
the homogeneous condition and the additional increase 
or decrease in sound pressure due to the inhomogen- 
eous condition. Using this concept, a meteorological 
function is defined which accounts for the departure of 
the atmospheric medium from homogeneous conditions 
and is given as 

>P = 7 (2) 

* 

where I* is the sound intensity received under homo- 
geneous conditions (obtained from equation 1 by noting 
that L;c = S*/pc). I is the sound intensity received 
under inhomogeneous conditions. 

A similar approach is followed by Brekhovskikh 
[3] in which ip was initially defined as a "focusing 
factor" and restricted to values of ip > 1. In general, 
ip can take on all values. Values of ip in the range 
0 ^ ip < 1 indicate a condition which is analogous to 
a "divergence zone"; i.e. , mean squared sound pres- 
sure is less than that under ideal conditions. Ideal 


in an inhomogeneous medium, the intensity I can be 
evaluated at any point along the tube. For the purpose 
of this analysis, the ray tube will be confined to lie in 
the plane which is defined by the y axis and the azimuth 
angle a. The shaded area at point R, (RFGJ) , is par- 
allel to the x - z plane and is separated by a vertical 

distance H from the x - z plane. The angle 6 is de- 

xt 

fined as the acute angle formed by the intersection of 
the perpendicular to the plane JKLR at R and the plane 
FGJR. 

The area of the wave front contained with in the ray 
tube at R is 

AB R = (RJ)(LR) 
where 

RJ s X A a 
s 

LR = RF sin 0 R 

LR = AX sin Q n . 

s R 
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The variable, X^, is the horizontal source -to-receiver 

distance of the ray tube and can be expressed as a 

function of the initial angle, 0 Q , of the emerging ray,, 

i.e.,X = f(0 ) # When this is done LR can written 
s o 


LR = 



0 

o 


A0 sin 0__ • 
o R 


The cross-sectional area 



o 


of the ray tube at R is then 

A0 sin 0 (3) 

o R 


where the absolute value signs restrict the equation to 

positive areas or positive energy flux through A EL.. 

R 

The power radiated by the source through the ini- 
tial surface area is 


AA = I AB 
o o 


where 


where <Mf,j3) is now evaluated at 0, with 0 being de- 

/ H \ ** ** 

fined as tan” \~^~) r is the straight line distance 

s 

OR from source to receiver. 


( 


> 



cos 0 r 2 
i o 



\4>(f,P) e 


>* 

o , 

\ 

X 


sin 0^ 

hi 

s 

\ 


®i 




o 


The above equation is a completely general ex- 
pression for ip and it is valid as long as the fundamen- 
tal assumption of ray acoustics is satisfied. If it is 
assumed, as in normal practice [5] , that the return 

angle of the ray, 0 , is approximately the same as the 
R 

emerging angle of tie ray, 0 , and that r » H, i.e. , 
r ~ X^then equation 7 reduces to 





/a e \ 

( c °t 0 o j < 

o 

(ax 

\ s'e 


°0 


ll. 


AB = (r cos 0 Aa)(r A0 ) 
o o o o o 


( 8 ) 


= r 2 A 0 A a cos 0 


and 


I o = A(f) <P(i 9 p) e j2 7 rrJ 


(I is the intensity of the source of radius r .) 
o o 


Combining equations 1, 2, and 8, while remember- 
iag that r « X g , yields 


S < f ’“’V X s ) = — 7 


pc A(f) <t>( t,p) e 


2*X s y(f,X s ) 


<p(f,a,e o ,x s ) 


(9) 


The power through A B q becomes 

AA = \(f) {a0 q A a cos ej /2 tt . • (4) 

By accounting for the molecular absorption losses , 
and assuming that the actual arc length of the ray tube 
is approximately equal to the horizontal source-to- 
receiver distance X and redistributing the power AA 
at R, the acoustic intensity I at point R becomes 

A(f)<P(UP) e cos 0 

AA_ o 

1 AB / 8X \ I 

2»X s T( ( ,X s ) (— ) I si "«R< 5 > 

0 

o 


where S(f, Of, 0 Q ,X g ) is the mean squared sound pres- 
sure that is received at a remote location from a point, 
directional, acoustic source radiating in an inhomog- 
eneous medium. For a nondirectional source, equation 
9 can be shown to be identical to the results presented 
in FI, 2]. 

In examining equation 9, it can be seen that the 
sound pressure is a function of the initial angle 0 q . 

This implies that it is the sound pressure associated 
with an individual ray. To arrive at the true sound 
pressure, all rays arriving at a common location 
should be combined. This will be shown in the follow- 
ing paragraphs. 


For homogeneous conditions, I jle ,from equation 1, 
at point R is 

A(f)<P(f£) e 

I* = — (6) 

2? r r 2 y (f,r) 


As stated previously 

X = f(0 ) . (10) 

s o 

This expression is normally evaluated by considering 
that the variation in the propagation \elocity of sound, 
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i.e. , the thermodynamic speed of sound plus the wind 
component, is a function of altitude only. Following 
along this same line, and using the hypothesis of a 
layered medium, atypical plot of X g *f( 0 Q ) for a given 
a?, is shown in Figure 3, The discontinuities in Fig- 
ure 3 are due to the common hypothesis of a layered 
medium [5, 6] with abrupt changes in the propagation 
characteristics between the layers. 



FIGURE 3. SOURCE TO RECEIVER DISTANCE 
VERSUS EMERGING ANGLE 

It can be seen from Figure 3, that for the indi- 
cated horizontal source-to-receiver distance, X , 

there will be five rays with initial angles of 0 , 0 ... 

°l <>2 

6 returning. In theory each ray has an associated 

mean squared sound pressure and therefore these rays 
have to be combined to arrive at the resultant mean 
squared sound pressure. If the pressures from each 
incoming ray are assumed to be uncorrelated, then the 
mean squared pressure associated with each ray, with 
each ray, rather than the instantaneous pressure, will 
add directly, and equation 2 can then be written as 

pc A(f) <P( Up) Q 

S(f,a,X ) = ~ | ip (f,ot,0 ,X )(11) 

S 2u X 2 y (f ,X ) i=l °i S 

s ' 9 s' 


where 

if>( f,a,0 ,X ) is given by equation 8 and the symbol 
o, s 

^ implies a summation of all^the rays arriving at a 

common landing distance, i.e, , receiver. The sum- 
mation inequation 11 is defined here, however, to in- 
clude only "direct" returned rays and not any "reflec- 
ted" returned rays. With the aid of equation 11, the 
expected mean squared pressure per unit bandwidth, 

P 2 / cps , received at any far -field location from abroad 
band, point, directional, acoustic source radiating into 
an inhomogeneous atmospheric medium can be com- 
puted. 


As noted previously, when the changes in the 
properties of the atmosphere within a wavelength of 
the energy under consideration can no longer be con- 
sidered as negligible, other means of arriving at the 
mean squared sound pressure have to be employed. 
This is readily apparent upon examining equation 8 


because 



increases without bond. This does 


not imply an infinite intensity, but rather that the fun- 
damental assumption of ray acoustics has been vio- 

/ ax g \ 

lated [2], The point at which \ — — j = o, i.e, , 

i x ) = °°> is defined as a focus and requires 

°i s 

a treatment similar to that presented by Brekhovskikh 
[3] to obtain an estimate of the pressure field. Like- 
wise, areas in which no ray can be returned is a fund- 
amental fault of the ray trace method [2] and does not 
imply zero pressure field. These areas are called 
shadow zones, and sound will be propagated into these 
areas because of diffraction effects [ 11,12] . 


SECTION IV. INHOMOGENEOUS SOUND PRESSURE 
EQUATION IN SOUND PRESSURE LEVEL FORM 


Sound pressures in acoustics are normally speci- 
fied in decibel form, i.e . , sound pressure level. The 
sound pressure level is defined 

i& \ 

SPL = 20 logic l ) dB (12) 

where ^/p^is the time root mean squared pressure and 
P ref iS * he r0ot mean squared reference pressure 

which is 0.00002 newton/ ( meter) 2 , Equation 12 can 
also be written as 


SPL = I0 1og 10 



dB 


(13) 


Dividing equation 1 1 by the reference pressure squared , 

i.e., (P 2 ) , and using the relationship in equation 13, 

ref 

equation 11 becomes the sound pressure spectrum 
level, or pressure spectrum level, (PSL), in deci- 
bels, 
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dB (14) 


dB 


(19) 


PSL(f* ft, X ) a 10 log,, 


pcA(f)^M) V , vl 

2 * X s V o y{f ' X j \ icl ^ <f * * °i ,X 7. 


Equation 13 can be expressed in its standard form 
as follows: 


PSL(f,a,X ) - 10 1 og 10 


\ P Vpc / 

O 


10 log,, <p(fJ3) 


- 10 log,o (2 7 TX z ) - 10 log,, y(f,XJ 

o S 


10 log , 0 


l V X .> 

i±l i 


dB ( 15 ) 


Using the units of newtons, meters, and seconds, 
P ref 

the ratio for standard sea level conditions is 

pc 

» 10~ 12 watts, which is the standard power reference 
when the SI system is used. The first term in equa- 
tion 15 then becomes the sound power spectrum level 
of the source, i.e., ASL (f) indB. 

From the ideas presented in Section II, the direc- 
tivity function, now called the directivity index, is 

/P*< UP) 0 \ 

I. (f,/3) * 10 1og 10 ( — ] dB (16) 

'h \ P*(f) j 

which is also equal to 

D.I. (f,P) 0 = SPL (f,/3) g - SPL (f) dB (17) 
h h 


where SPL (f) is called the space average sound pres- 
sure level. 


Similarly, from Section II, the attenuation func- 
tion, now called the attenuation factor, becomes 

E A(f,r) = 10 lo eto Y( f .r) = 10 log 10 e 2 * (f,r . (18) 

By further reductionand noting the relationship between 
nepers and decibels and also the fact that r » X , the 
above becomes s 


E A(f,X s ) 


r (f) X s 


where e (f) is now the attenuation coefficient in 
dB/meter and X is the horizontal source to receiver 
distance in meters. 


The last term in equation 15 is called the ’’mete- 
orological factor” and is 


M (f, o,X ) 4nT 
F s' = 10 Iog 10 


t, -V 

1=1 


dB. (20) 


Equation 15 now becomes 


PSL(f,Q,X ) = ASL(f) + D, I. (f,/3) -20 1og 10 X 

^ s 

n 

- E A < f > x s ) + M F (f,a,X s ) -8 dB. (21) 


Due to the broad band frequency characteristics 
of the source, pressure spectrums in sound pressure 
level form normally have associated bandwidth charac- 
teristics. 

In the case of far -field acoustics, the bandwidth 
associated with the pressure spectrums are usually 
octave bandwidths. Equation 21 is related to an octave 
bandwidth analysis as follows 

PSL(f,a,X s ) + l01og 10 Af=O. B. SPL(f,Of t X ) dB 

( 22 ) 

where Af = 0.707 f , with f being the center frequency 
c c 

of the octave band of interest. Equation 21 then be- 
comes 

O. B. SPL(f,a,XJ = O. B. AL(f) +D.I.(f ft) 

S °h 

-20 logjo X g - E A (f,X g ) + M p (f, Q,X g ) - 8 dB(23) 


where 

O. B. AL(f) = ASL(f) + 10 log 10 Af dB . 

Equation 23 is then the standard form in which t’ 
far-field sound pressure level prediction equation is 
normally presented. 



SECTION V. DISCUSSION OF THE FAR -FIE LD 
SOUND PRESSURE LEVEL EQUATION 


In "rocket engine generated acoustics, " the stand- 
ard or conventional form of the equation that is used to 
obtain an estimate of the sound pressure level, in deci- 
bel (dB) , at a remote far-field location resulting from 
the operation of a rocket engine, is that which is pre- 
sented by equation 23. Combining equations 8, 20, and 
23 this equation is again presented here and is 

O. B. SPL(f,a,X )= O. B. AL(f)+D.I. (f,/3) 

8 h 

-20 1og 10 X g -E A (f,X a ) -8 


meteorological factor of equation 24 is made to b 
unity, this would make the independent of any 
source characteristics and would allow to be eval- 
uated solely in terms of the existing atmospheric con- 
ditions. The amount of error introduced in the com- 
puted O. B. SPL, by making the directivity ratio unity, 
is not presently known. This idea will be investigated 
and if the errors are small, then it would seem logical 
to evaluate for many different recurring atmos- 
pheric conditions [13]. 

In any event, equation 24 can be used in its exist- 
ing form to obtain engineering estimates of the sound 
pressure level received at a remote far-field location 
from a broad band, point, directional, acoustic source 
radiating into an inhomogeneous medium. 


o. B. SPL<f,<*,X g ) = 0. B. AL<f) ♦ D.I. - 20 log], X g - E^f.X^ - 8 


■ 10 log], 


I ( X -K») s l( -, 4 «'■<)} 


In examining equation 24, it can be seen that the 
first five terms on the right side of the equation can be 
combined to yield the octave band sound pressure level 
that is observed under the ideal conditions of a homo- 
geneous atmosphere, i.e. O. B. SPL (f, a, X g ) . 

Equation 24 can therefore be written as 

O. B. SPL(f,o,X ) =0. B. SPL(f,a,X ) 
s s J e 


+ M (f, o,X ) dB, 

r S 


(25) 


From the above, it is obvious that the "meteoro- 
logical factor," i.e. , M (f,a,X ), describes the ef- 

r S 

fects of the inhomogeneity of the atmosphere on the 

received O. B, SPL. This term in essence "calibrates” 

the atmosphere. The unfortunate part, however, is 

that the IVT term is a function of the directional char - 
F 

acteristics of the source under consideration. 

If a nondirectional source is considered, the me- 
teorological factor becomes independent of the source 
characteristics (eq. 24) , and is in turn solely a func- 
tion of the sound propagation velocity. For this case, 
the meteorological factor could be evaluated for many 
different types of atmospheric conditions without having 
any prior knowledge of the characteristics of the 
source. 

A nondirectional source rarely ever exists in ac- 
tual practice. However, if the directivity ratio in the 


SECTION VI. CONCLUSIONS 


From the ideas concerning acoustic source char- 
acteristics and ray acoustics discussed thus far, an 
expression has been derived, i.e. , equation 24, which 
is of extreme importance is the area of "far-field 
acoustics.” This equation has a direct application in 
far-field sound pressure level studies in connection 
with the static firing and launch of any space vehicle 
system. Equation 24, although very useful, is by no 
means the final solution. As stated previously, the 
’’meteorological factor” has a severe limitation which 
results from the failure of ray acoustics, in that no 
finite solution can be obtained at a focus point. An 
alternate approach has to be used [3] Li order toobtain 
a solution in these regions. This factor, i.e., M^, 
appears to be extremely sensitive in only the imme- 
diate area of a focus (Fig. 3) . Another , but not as 
severe, limitation of the term, is that no solution 
can be obtained in the regions where no rays are re- 
turned, i.e. , shadow zones. Even though these regions 
indicate that no sound rays are present, sound will pro- 
pagate into these areas because of scattering and dif- 
fraction effects [11, 12]. 


In order to use the far-tield sound pressure level 
equation, three quantities have to be known. The first 
two, i.e., ORAL (f) and D.I. (f,/?) , are concerned 
with the acoustic source in question. These quantities 
are normally determined by model tests or some re- 
liable empirical method. The third quantity (X = 
{(0 )) is somewhat more difficult to evaluate. As in- 
dicated earlier, the normal approach, and by far the 
simplest, is to assume that the sound velocity can be 
represented by a single profile, which is a function of 
altitude only, and further, that the linear layer hypo- 
thesis is valid. The validity of both assumptions is 



questionable especially the latter one. It has been 
indicatecrby several authors [14] that very severe 
anomalies result because of the use of a linear layer 
model. To mention a few, focal zones appear where 
they should not, and where these zones should be pres- 
ent, none is indicated. Numerous studies involving 
very elaborate curve-fitting routines have been per- 
formed [15, 16] in an attempt to alleviate this as- 
sumption, The necessity in the introduction of many 
different vertical velocity profiles as a function of 
horizontal distance [4] would depend upon local en- 
vironmental conditions.. In some areas, the normal 
assumption is probably adequate. Equation 24 imposes 

no restriction on how (X = f( 6 ) ) is to be evaluated. 

s o 

The only restriction is that the sound rays, for a given 
azimuth, is confined in a vertical plane. 

In summary, and with due recognition to the lim- 
itations noted above, it is believed that reliable en- 
gineering estimates of the sound pressure received at 
remote far-field locations from a broad band, direc- 
tional, point source radiating into an inhomogeneous 
medium can be made with the aid of equation 24 and a 

functional form of X = f(0 ), similar to that which 
s o 

is presented in reference 5. 
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LAUNCH VEHICLE NOSE SHROUD OPTIMIZATION 
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SUMMARY 

2 - i 

The effect of nose fairing shape on vehicle per- 
formance has been studied for two -stage and three- 
stage Saturn V type vehicles. Shroud shapes con- 
sidered were single and double angle cones with and 
without cylindrical extensions. Other factors that in- 
fluence shroud weight and vehicle performance were 
also investiga ted. The more important variables were 
found to be the shroud shape, time of jettisoning, re <■ 
quired volume, choice of material, and type of con- 
struction. The trends found are generally applicable 
for any vertically launched space vehicle, having a 
thrust-to-weight ratio at lift -off of approximately 1, 25. 

'0 >jT : 

I. INTRODUCTION 


A parametric study has been conducted to deter- 
mine the various factors that influence the design of 
nose shrouds for large launch vehicles and todetermine 
the effect of these factors on vehicle performance. 2 
The study was initiated by MSFC early in 1962 in re- 
sponse to a request from the Jet Propulsion Laboratory 
to determine from a performance standpoint the most 
profitable means to package payloads in Saturn -class 
vehicles. The results of the study can also aid in the 
selection of the most appropriate nose shape to shroud 
a given payload. 

Fairing shapes considered were single and double 
ingle cones, with and without cylindrical extensions. 


1 Research Specialist, Lockheed Missiles and Space 
Company, Huntsville Research and Engineering Cen- 
ter, Huntsville, Alabama, 

2 The study summarized in this report was performed 
by LMSC/HREC for NASA/ MSFC under Contract 
NAS8-11148. 


Although the majority of the study was conducted using 
a three -stage Saturn V vehicle (Fig, 1), payload 



weights were also determined for 
a two-stage Saturn V vehicle to de- 
termine the dependence o. the ”op- 
timum" nose shapes on the launch 
vehicle and mission selected. Oth- 
er r ictors considered in tile study 
wei- shroud material and type of 
construction, nose bluntness, jet- 
tison time, aerodynamic drag, 
aer ody nam ic 1 oad s , and ae r ody - 
namic heating. 

A typical double angle nose 
shroud showing the basic dimen- 
sions iS shown in Figure 2, Initial 
cone angles (flj) used were 15°, 
20°, 25°, and 30°; frustum angles 
(0 2 ) were 5°, 7.5°, 10°, 15°, and 
20°; and diameter ratios (D { /D 2 ) 
were 0.4, 0.6, and 0.8. In addi- 
tion, single angle cones having half 
ahgles from 10° to 25° were studied. 


FIGURE 1. S \TURN V THREE ST ,GE VEHICLE 
WITH TYPICAL NOSE SHROUD 



FIGURE 2. SCHEMATIC OF TYPICAL NOSE SHROUD 
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II. METHOD OF ANALYSIS 


The performance of a launch vehicle is both di- 
rectly and indirectly a function of the shape of the nose 
fairing which, in most cases, also serves as a pro- 
tective shroud for some type of payload. Performance 
is directly affected by the influence of aerodynamic 
drag which is dependent on the nose shape. Indirectly 
the performance is affected in that the shape influences 
parameters such as aerodynamic loads and heating, 
which in turn influence the weight of the structure. 

Numerous investigations have been made to de- 
termine minimum drag shapes using the calculus of 
variations in conjunction with various theoretical 
means. In general, for a body of given fineness ratio, 
the minimum drag shape is described by some ’’power 
law. " However, fabricating a contoured shroud of this 
type would be highly impractical for the sizes required 
for large launch vehicles. In addition, the minimum 
drag configuration does not necessarily give the most 
performance due to structural weight considerations. 

Assuming, then, that an "optimum” shroud would 
require ease of fabrication, body contours would have 
to be approximated by a number of straight line seg- 
ments. To keep the number of variables to a minimum, 
initial studies considered only conical shapes fl] and 
were later extended to double angle ( cone plus frustum) 
configurations [2], 

By choosing various cone angles, frustum angles, 
and diameter ratios, the aerodynamic drag, loads, 
and heating were determined [3, 4, 5]. The structure 
of each shroud was then designed using computer pro- 
grams, developed during the course of the study, which 
minimize shroud weight [6,7] . 

Fitting the vehicles with various nose fairings and 
holding stage gross weights constant, ascent trajec- 
tories into a 185. 2 kilometer circular orbit were op- 
timized. Afterwards, the three -stage vehicle con- 
tinued on an earth-escape mission. The effect of 
jettisoning the nose fairing at different times during 
flight was co* r r red for both vehicles and missions. 
From the resulting trajectories, maximum payload 
capabilities were determined. 

in. RESULTS 


Of the parameters studied, the me having the most 
influence on performance is the shroud jettison time. 
For the Saturn V three -stage -to -escape vehicle, the 


payload can vary as much as 5,000 pounds {2270 kg) 
(approximately 5 percent) depending on whether the 
shroud is jettisoned after first stage flight or is car- 
ried toeseape conditions. Not only does jettison time 
influence payload weight, but it also greatly influences 
the "optimum" shroud shape. Jonsidering a conical 
nose shroud, it can be seen that this is primarily a 
result of the opposing influence of cone angle on drag 
and fairing structural weight. As shown in Figure 3, 



FIGURE 3. EFFECT OF CONE ANGLE ON TOTAL 
POWER-ON DRAG COEFFICIENT 


the total vehicle drag increases with increasing cone 
angle. However, in Figure 4, it is shown that as the 
cone angle increases, the shroud structural weight 
decreases until the cone angle reaches approximately 
25°. As would be expected, the resulting performance 
curves for conical shrouds (Fig. 5) show that when 
the shroud is ejected early in flight, the aerodynamic 
drag is of primary importance and a shallow angle 
shroud gives the most performance. The longer the 
shroud is carried, the more important weight be- 
comes, and thus the maximum performance cone angle 
is somewhat higher. 

For a given ejection time, varying the shroud 
shape gives differences in maximum payload of from 
1,000 to 4,000 lbm (454 kg to 1814 kg), indicating that 
nose shape is also an important factor in vehicle per- 
formance. Figure 6 shows typical results for double 
angle fairings having an initial cone angle of 25° and 
diameter ratio of 0.4. For early jettison times, 
shallow angle frustums appear to give the most 
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EFFECT OF CONE ANGLE ON TOTAL 
CONE WEIGHT 
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FIGURE 6. EFFECT OF FRUSTUM ANGLE ON 
PAYLOAD WEIGHT 



performance, whereas when the shroud is carried to 
escape, the most performance is obtained when the 
frustum angle is closer to the initial cone angle . Again 
note the marked effect of jettison time on payload 
weight. 

It is to be expected that many low density payloads 
cannot be contained within the volume of the double angle 
cone giving the most performance. Thus, to determine 
the best means of acquiring additional volume, the ef- 
fect of cylindrical shroud extensions on payload weight 
was established. A typical result from Reference 2 
is shown in Figure 7 where maximum payload weight 
is presented as a function of the shroud gross volume 
and geometry. These curves indicate that, for the 
conditions given, it would be more advantageous not to 
use cylindrical extensions if the required volume is 
less than approximately 10, 000ft 3 (283 m a ). However, 
for greater volume requirements, the use of cylindri- 
cal extensions would give the greatest payload capa- 
bility. 


Other factors, mainly, shroud material, type of 
FIGURE 5, EFFECT OF CONE ANGLE ON construction, and internal pressure variations from 

VEHICLE PAYLOAD COMPABILITY ambient, were ound to have a pronounced effect on 
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shroud weight. However, to keep the number of vari- 
ables within reasonable limits for the time allotted for 
the study, the bulk of the present w or k considered only 
aluminum ring plus skin designs and ambient internal 
pressure. 



Notes: I. Three -Stage SATURN Z Vehicle 
90 [~ 2. Fairing Jettisoned after Second Stage Flight 

3. Aluminum Shrouds 


89 Y 0| = 25* 

D,/0 2 * 0.4 
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FIGURE 7. EFFECT OF FRUSTUM ANGLE AND 
TOTAL FAIRING VOLUME ON PAYLOAD WEIGHT 

An enlightening result of the studies thus far is 
the fact that the same trends of payload weight versus 
body geometry are evident for both the two -stage and 
three-stage vehicles investigated regardless of their 
missions. It can thus be conjectured that the trends 
observed are generally appropriate for any vertically 
launched space vehicle having an initial thrust-to- 
weight ratio of approximately 1,25 and similar dy- 
namic pressure histories. 


IV. CONCLUDING REMARKS 


The nose shroud optimization studies conducted 
thus far have determined: 

a. The factors that most influence shroud weight 
and vehicle performance. 

b. Within limits, the most appropriate payload 
shroud shapes for maximum performance, 


c . The areas where additional effort is necessary 
and/or desirable. 

The factors having the greatest effect on perform- 
ance were found to be shroud shape, time of jettison- 
ing, required volume, choice of material, and type of 
construction. Although the two-stage and three -stage 
Saturn V vehicles gave considerably different payload 
weights, the same trends have been observed for both 
vehicles. Theresultscanthusbeconsidered generally 
applicable for any vertically launched space vehicle 
having a lift-off acceleration of approximately 1.25 and 
similar dynamic pressure histories. 

In most instances, especially when the shroud is 
jettisoned early in flight, more performance can be 
obtained by the use of double angle cones as opposed 
to pure conical shrouds. In addition, the double angle 
cones have more volume per unit length for better 
payload packaging. 

Although a great deal of useful information has 
been derived from the studies thus far, there are sev- 
eral areas where additional effort should prove par- 
ticularly rewarding. One of these, structural design 
optimization, is being explored by the Lockheed Mis- 
siles and Space Company, Huntsville Research and 
Engineering Center for the Propulsion and Vehicle 
Engineering Laboratory at MSFC. The present work 
has been based upon weights of optimized (weight-wise , 
subject to practical constraints) ring plus skin de- 
signs. Two alternate methods of fairing construction 
that may be more desirable for certain applications 
are ring plus stringer plus skin and honeycomb sand- 
wich. Incorporating these materials into a shape 
optimization study will require that similar structural 
optimization computer programs be developed. 

Since it has been found that the internal pressure 
in a nose fairing can significantly affect the design of 
the fairing, an in-flight venting analysis should be per- 
formed and included in the design optimization studies. 
Computer programs are presently available that pre- 
dict the variation of internal pressure with flight time 
provided the external flow properties at the vent loca- 
tion are known. 

Another area where additional information is 
needed is in pressure and load distributions at trans- 
onic Mach numbers. The present studies have been 
performed using a set of parametric charts [ 4] derived 
from a very limited amount of test data. With addi- 
tional experimental results , more reliable charts could 
be developed tAus enabling closer approximations to 
the true minimum weight shroud designs. 
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There is also a possibility that shapes other than 
single and double angle cones would be more desirable. 
The payload gains from parabolic or other contoured 
configurations could conceivably offset the difficulties 
of fabrication and thus should be considered in future 
studies. 
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TRISONIC AERODYNAMIC LOADS STUDY 

By 

Bob G. Dunn 


SUMMARY 

This study is intended to define the linear aerody- 

Symbol 

Definition 

M 

Pitching moment 

namic load distributions on basic missile and launch 
vehicle shapes. Since theoretical means are not avail- 

M c 

Critical Mach number 

able in the transonic regime to accurately predict loads , 
experimental data are being used in the Mach number 

M 

x> 

Freestream Mach number 

range from 0,7 to 2, 0. 

p 

Surface static pressure 



yo 

‘ 


By correlating data from the literature with the 
results of a supplementary test program, aerodynamic 
load distributions for cones and cone -cylinders have 
been obtained and published. 


Free stream static pressure 


Freestream dynamic pressure, 
0.7 p^AI^ 2 


A test program is in progress to define linear 
local loads on cone -cylinder -flare -cylinder shapes. 
The analysis of these data has recently begun, the final 
results of which are to be published in August I960. 

- - 


s 

O' 

o 


Reference area, ?r D 2 /4 
Angle -of-attack, degrees 
Cone angle, degrees 


LIST OF SYMBOLS 


Subscripts 


Symbol Definition 

C^ T Normal force coefficient, F XT /q S 

N N 

C >T Normal force coefficient slope, 

N « 0 C N /0 « 

C ' Local normal force coefficient, 

9 C N /0 (x/D) 

C ' Local normal force coefficient sloix', 

N « 9 C^/9 cv 

C Pitching moment coefficient, 

M M/q SD 

C Pitching moment coefficient slope, 

« 9 C m /9 o- 

C Pressure coefficient, (p-p^)/q^ 

D Cylinder diameter 

F XT Normal force 

N 

L Segment length (n denotes segment 

number) 


1 Cone segment 

2 Fore -cylinder segment 


I. INTRODUCTION 


Missiles and space vehicles are often configura- 
tions made up of several conical and cylindrical seg- 
ments or variations of these shapes. Some examples 
of these are the Saturn vehicles, Minute man, Red- 
stone, Jupiter, Juno, Titan, Polaris, and Agena. 
For the trajectories flown by most of these vehicles, 
maximum or "design" aerodynamic loading conditions 
occur at speeds near Mach one. 

Even though much testing and analysis are done 
after initial design to determine exact loading and per- 
formance characteristics, the vehicle design and fab- 
rication are commonly ha>ed upon aerodynamic loads 
and characteristics established in the preliminary de- 
sign analysis. Thus, the importance of having good 
preliminary design cannot be over -emphasized. 

Although theoretical methods are quite adequate 
for preliminary design in other areas, transonic flow 


17 



fields about these configurations do not lend them- 
selves to accurate theoretical analyses. This being 
the case , the designer must rely upon experimental 
aerodynamic design curves for initial design infor- 
mation. A parametric representation of these factors 
is a prerequisite to any practical optimization studies. 
Such a practical set of design curves which cover 
broad variations of the significant parameters have not 
previously been available. 

For the past several months, the Huntsville Re- 
search & Engineering Center of Lockheed Missiles & 
Space Company has been engaged in a program to gen- 
erate an adequate set of aerodynamic design curves for 
basic shapes in the "trisonic" Mach number regime 
(high subsonic, transonic and low supersonic). The 
Trisonic Aerodynamic Loads Study is under the direc- 
tion of the Aero-Astrodynamics Laboratory of Marshall * 
Space Flight Center, National Aeronautics and Space 
Administration. 


II. DISCUSSION 

A, Published Results 

A "building-block'' approach to this task was 
adopted whereby the simplest configuration, a cone, 
was considered first, and more complicated bodies 
were considered in the order of complexity of shape. 
It w as planned that the effects of each parameter on the 
loads for a given segment would be defined such that 
the total loads for a composite body could be obtained 
by a summation of segment loads. 

A search of the literature was performed to collect 
published information pertinent to the study. Early in 
the search, it became apparent that the quantity of use- 
ful experimental results was limited and that a test 
program would be necessary to provide the data re- 
quired for a comprehensive study. The correlation of 
available results for cones, cylinders and flares, pub- 
lished in Reference 1 , did not provide continuous rep- 
resentation of pressure distributions and stability 
characteristics for the area of interest. 

To define the aerodynamic loads on cone -cylinder 
configurations, a series of wind tunnel tests was con- 
ducted to supplement existing results, thus providing 
a continuous spectrum of data for the appropriate range 
of each parameter. Both force and pressure tests 
were conducted in the MSFC 14 x 14 inch Trisonic Wind 
Tunnel. A subsequent analysis combined collected 
results from the literature with those of the associated 
test and resulted in the design curves presented in 
Reference 2. The linear aerodynamic characteristics 
of cones and cylinders are shown in the form of 


(1) local normal force coefficient slope, C< T 

iN o> 

(2) surface pressure coefficient, C p 

(3) normal force coefficient slope buildup, C*-, 

1N 0! 

(4) pitching moment coefficient slope buildup, 



The parameters and ranges of each are 

(1) Mach number (M m ) from 0.7 to 2.0 

(2) cone semi -vertex angle (0j) from 10 to 40 
degrees 

(3) cylinder lengths up to 6.0 diameters. 

Some typical curves are shown in the "carpet 
plots" of Figures 1 through 4, the experimental points 
shown being taken from References 3 through 6. The 
carpet plot method presented here was chosen for two 
reasons: (1) the coefficients presented are shown 

varying continuously with two parameters, and (2) 
for correlation purposes , points which do not fall on 
one particular curve may be used to define the whole 
"family" of curves as opposed to the requirement that 
all parameters but one be duplicated. 

Figure 1 presents the variation of cone surface 
pressure distribution with axial station and cone angle 
for M^ =0.8. A complete cone C p distribution for 

any cone angle between 10 and 40 degrees can be faired 
in by connecting the proper points on each constant sta- 
tion curve. 

The cylinder surface pressure for cone -cylinders 
is shown as a function of axial station and cone angle 
at Mach 0. 8 in Figure 2. It is seen that there is a 
large pressure drop due to the expansion around the 
cone -cylinder junction, particularly for the smaller 
cone angles. The sharp break in the x/D = 0. 0 curve 
just above 20° is characteristic of boundary lay«r sep- 
aration at the juncture. For larger cone angles, the 
flow is separated, whereas it remains attached for 
smaller cone angles. The Mach number correspond- 
ing to the cone angle whore separation occurs is term- 
ed the critical Mach number and will be discussed 
later. 

Figure 3 depicts the local normal force slopes 
for cones as being dependent on cone angle and axial 
station. Integrated values for these curves are com- 
pared with cone force data in Refex'ence 2 and are 
in excellent agreement. It is signiiicant that the pre- 
sence of a cylindrical afterbody makes little difference 
on the cone loading. 
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FIGURE 3. THE EFFECT OF SEMI-VERTEX ANGLE ON LOCAL NORMAL 
FORCE SLOPE DISTRIBUTIONS FOR CONES 

The local normal force slopes for the cylinder seen in the C p curves. A more detailed survey of data 

portions of cone-cylinders are shown in Figure 4. In in the vicinity of 0= 20° might indicate some abrupt 

this case, one sharp break appears in the constant sta- change in the trends as was observed in some other 

tion curves (at 6 = 20°) due to flow separation as was cases. 



FIGURE 4. THE EFFECT OF CONE ANGLE ON LOCAL NORMAL FORCE 
SLOPE DISTRIBUTIONS FOR CONE -CYLINDERS 
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T 1 ' curves for each segment may be combined to 
obtain ,ne load distributions for a cone -cylinder; some 
typical cases are presented in Figures 5 and 0. The 
pressure distributions are shown in Figure 5 and local 
normal force slopes in Figure 6 for Mach numbers 0.8, 
1. 1, and 1.5, The curves for Mach 0.8 are taken di- 
rectly from interpolated curves faired in Figures 1 
through 4, These Mach numbers represent three dis- 
tinct flow patterns : Meo = 0,8 is a transonic separated 
tlow case, M„ = 1. 1 is a transonic case with boundary 
layer completely attached, and Mco = 1.5 is a super- 
sonic case. 



FIGURE 5. SURFACE PRESSURE DISTRIBUTIONS 
FOR A 22,5 DEGREE CONE-CYLINDER CONFIG- 
URATION AT ZERO ANGLE-OF-ATTACK AND 
MACH NUMBERS 0.8, 1.1 AND 1.5 

For zero angle -of -attack and fixed cone angle, 
there is an approximate critical Mach number, M^, 
above which flow is attached and below which the 
boundary layer separates at the juncture. Very close 
to this critical Mach number, alternating flow attach- 
ment and separation may occur, but such a phenomenon 
has not been explored in the present study. Reynolds 
number and possibly other parameters may have con- 
sklerable effect upon the critical Mach l umber. How- 
ever, this is a problem worthy of a separate research 
and is not investigated here. 

To define the critical Mach number well, the 
models should be tested at many Mach numbers near 




Model Station* X/D, callbere 


FIGURE 6. LOCAL NORMAL FORCE SLOPE DIS- 
TRIBUTIONS FOR A 22.5 DEGREE CONE- 
CYLINDER CONFIGURATION AT MACH 
NUMBERS 0.8, 1.1 AND 1,3 

Mq, which has not been done for the present investi- 
gation. It was observed from the data available, how- 
ever, that both separated and attached boundary layers 
existed for given conditions. The Mach numbers for 
such cases were inferred to be Mq. Figure 7 shows 
the variation of M^ with cone angle as derived in this 
manner compared with similar results acquired from 
the pressure curves of Reference 4 . The curves do 
not agree, but do converge for higher Mach numbers. 
The discrepancy is probably due to the different wind 
tunnels or other unexplored factors. 

B. Current Activities 

The Trisonic Aerodynamic Loads Study is be- 
ing extended to include the effects of flare -cylinders 
attached downstream to cone -cylinders. Since the 
literature provided very few test data which are ap- 
plicable, a complete test program was scheduled in 
the MSFC 14-in. (3(> cm) Trisonic Wind Tunnel toob- 
fcain pressure distributions on cone -cylinder -flare- 
cylinder configurations. A force was also scheduled to 
verify integrated pressure loads. Details of these test 
plans are documented in Reference 7 , 


^1 



FIGURE 7. THE INFLUENCE OF CONE ANGLE ON 
CRITICAL FLOW SEPARATION MACH NUMBER 

At the present time , a large portion of the pres- 
sure test has been completed, A second entry is 
scheduled for September 1965 to complete the pro- 
gram; the force testwill also be completed at that time . 

Analysis of tb"! test results generated in the test 
mentioned above ha3 recently begun. It is the purpose 
of the current phase of study to derive and present 
parametric design curves which define linear local 
aerodynamic load distributions for cone -cylinder - 
frustum -cylinder bodies, interest will center about 
the frustum, and an attempt will be made to isolate 
the flare effects from the nose shape effects. In this 
manner, the building -block approach can be followed 
in the use of these and previously published curves 
provided judicious restraint is practiced in their ap- 
plication. 

The Mach number range covered will be from 0.7 
to 2.0, Flare angles of 7.5°, 10°, 15° and 20° and 
flare base -to-leading -edge diameter ratios of 0.4, 

0.6, and 0.8 will be considered in the coming analysis. 

It is expected that the results of the current phase 
of 3tudy will be published in a handbook by August of 
1966, The planned curves will provide the preliminary 
designer with a convenient means of evaluating the 
aerodynamic loads for configurations falling into this 
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general category without time-consuming and expen- 
sive tests. The optimization of various parameters 
with respect to vehicle performance will also be prac- 
tical with the availability of the final results. 
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ON SOME AERODYNAMIC ASPECTS OF THE PANEL FLUTTER PROBLEM 



By 

M. PJatzer, R. Beranek, L. Saunders 


<r . 



SUMMARY 


Symbol 


Definition 


Recent panel flutter studies have led to the con- 
clusion that linearized potential flow theory is inade- 
quate to predict panel flutter at low supersonic Mach 
numbers. Nonlinear transonic flow effects and bound- 
ary layer effects were suggested as major reasons for 
the disagreement between theory and experiment. This 
paper presents a comparison between linear and non- 
linear invisc id supersonic flow past wavy walled planar 
and cylindrical configurations and discusses their 
significance toward panel flutter research. 

LIST OF SYMBOLS 



s 

t 


u, v,w 


x,r ,0 

a 


a" 


Entropy 
Ti. . 

Velocity components in the x, r , 6 di- 
rections of the velocity V 

Coordinate system 

Wavelength of the advancing acoustic 

disturbance 777777 
c (M+l) 

Wavelength of the^receding acoustic 
disturbance 


c Velocity of sound in the free stream 

0 ^ Specific heat at constant pressure 

c v Specific heat at constant volume 

C Euler *s constant = 0.5772157 

J (A) Bessel function of the first kind of order n 
n 


y Wave number ( 27 r/L) 

£ Real , positive variable 

c p Velocity potential 

co Circular frequency of the cylinder 

vibration 

v Ratio of specific heats 


Y (A) Bessel function of the second kind of order n 
n 

( 2 ) 

H q (z)Hankel function (Bessel function of the 
third kind) 

K (z) Modified Bessel function of the second kind 
o 

L Wave length 

M Mach number 

n Direction normal to a streamline 

R Gas constant 

R Radius of cylinder 

o 

R Surface deformation of the cylinder, posi- 

tive in the positive r direction 

R(x,t) Radial displacement of the cylinder 


I. INTRODUCTION 


The dynamic aeroelastic instability of surface 
skin panels that have one side exposed to an air stream 
and the other side to still air has become of increas- 
ing concern to designers of air craft and missile struc- 
tures in the last ten years. Tins problem is commonly 
referred to as panel flutter and is defined as the self- 
excited oscillation of a panel in a flow under the action 
of aerodynamic forces generated by the panel motion. 
Experiments have shown conclusively that flutter of 
surface skin panels can exist [ 1, 2 ], As a matter of 
fact, Jordan [3] has suggested that such flutter con- 
tributed to early German V -2 rocket failures. More 
recently, panel flutter has occurred on the X-15 air- 
plane. Figure i shows those regions on the X-15 
re sear chair craft which were affected by panel flutter. 
The variation of the panel response envelope with dy- 
namic pressure as obtained from flight measurements 
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devoted to this problem in the last decade, the areas 
of agreement between theory and experiment are quite 
limited. As a result of these studies (see Ref. 7 for 
a recent review) it is generally realized now that in 
panel flutter many factors are important, e. g. , the flow 
conditions, initial curvature effects, tension and com- 
pression loads induced by mechanical or thermal ef- 
fects, etc. Thus, a comparison between theory and 
experiment may allow proper interpretation only if all 
these parameters have been thoroughly analyzed and 
controlled, an obviously difficult task. 


FIGURE 1. REGIONS OF X-15 RESEARCH AIR- 
CRAFT AFFECTED BY PANEL FLUTTER [ 5] 

on an X-15 side-fairing panel is shown in Figure 2. 
The upper curve exhibits a sudden increase in ampli- 
tude at a dynamic pressure of 600 psf due to start of 
flutter. The lower curve reproduces the measured 
envelope of the response for the stiffened panel, which 


Further progress in the analysis of panel flutter 
therefore seems to depend upon a careful separation 
and evaluation of the different factors influencing panel 
flutter experiments. In this respect, Lock and Fung 
[81 have suggested that the effects of aerodynamic 
nonlinearity and of viscosity may be quite important. 
These aerodynamic aspects of the panel flutter problem 
are discussed in this report in more detail, along with 
some results of preliminary in-house studies. 


RELATIVE 

PANEL 

RESPONSE 



FIGURE 2. VARIATION OF PANEL-RESPONSE EN- 
VELOPE WITH DYNAMIC PRESSURE OBTAINED 
FROM FLIGHT MEASUREMENTS ON AN X-15 
SIDE-FAIRING PANEL [4] 

shows that flutter was eliminated. Although, in gen- 
eral, structural nonlinearities strongly tend to limit 
the flutter amplitudes and therefore cause the modes of 
structural failure to be those peculiar to fatigue rather 
than explosive fracture of the skin surface, it is dan- 
gerous to regard panel flutter as a fatigue problem 
only. Tests of the corrugation stiffened X-15 panels, 
e. g. , revealed that flutter was very sudden and very 
severe. The panels failed unless flutter was stopped 
in a few seconds [4] . Similar explosive type flutter 
was encountered in tests of corrugation-stiffened 
Dyna-Soar panels [ 6] . 


II. PRESSURE DISTRIBUTION ON A VIBRATING 
CYLINDER 


We select a vibrating cylinder for a discussion of 
the physical significance of various aerodynamic the- 
ories. Consider an infinitely long cylindrical shell of 
radius R 0 which is exposed externally to a uniform air- 
stream parallel to the generators of the cylinder. As- 
sume the cylinder to be vibrating in a sinusoidal axi- 
symmetric standing wave mode within the length L. 
Thus, the panel deflection is given by 

R (x, t) = R sin (y x) • , ( 1) 


where R (x, t) is the radial displacement of the shell 
measured positive in the outward direction, R is a 
real or complex number, and y is the wave number of 
the standing wave deflection. The deformation of the 
shell is assumed to be so small that’ linearized aero- 
dynamic theory may be employed. 


Therefore, adopting a cylindrical coordinate sys- 
tem, we are looking for a solution of the linearized 
unsteady potential equation. 


...» ... 1 . . 2io>M , o> 2 , A 

(M 2 - i) <p - $ “ 7 </> + “ 0 (2) 


xx 


rr 


x c‘ 


Unfortunately, there is no reliable formula or 
computing program of guaranteed accuracy available 
today to guide missile and aircraft engineers in de- 
signing against panel flutter. Although a substantial 
amount of theoretical and experimental effort has been 


Such a solution may be obtained by distributing 
pulsating sources along the axis and summing the ef- 
fect of the sources. The basic source solution in 
supersonic flow is well known [ 9] ; e. g. , 
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1 


0 .(X, r,t) =- 


,i[ cot - M (*-£)] 


2 ir nT ( x-|) 2 - cot 2 a • r 2 
cos (W (x-£) 2 - cot 2 a • r 2 )', 


(3) 


to x - r cot a instead of integrating from the leading 
edge x = 0. This results in a considerable simplifi- 
cation because the integration over the x - direction 
can now be carried out and we obtain [ 11] 


where 


cot 2 a = M 2 - 1, k = 


co 


c cot 2 a 


i o M 
c cot 2 oi 


Thus, by summing over all these sources from 
0 to X, we obtain 




, i icot -iyx r _ . 4 TT ( 2 ) 

+ 4 e q 2 e {P 2 (r) - IT 1 


where 


< r ) “ Hq 2) (fr)} 
(~r)} (7) 


<P(k, r,t) = 


-icot x-r cot a f , $ 5 5-1] 

\ r q(£) cos |/W (x-S) -cot o? • r_ ) 

0 \f(x-|) 2 - cot 2 a • r 2 * 


Pj(r) = [ r cot crJ( 0 !+y)((J!" + y)'] (7a) 

P 2 ( r ) =H Q 2) tr cot a 4(5 - y) (a" - y)'] ^ 




d{ • 


(4) 


The source distribution can be determined from 
the boundary conditions at the cylinder surface, and 
in our following considerations, we regard it to be 
known. 


Equation (4), unfortunately, is unsuitable for a 
physical interpretation of the generated flow field. 
We, therefore, try to generalize Munk f s cross-flow 
concept [10] to this case. Thus, considering first 
only the flow in a cross section perpendicular to the 
free-stream, we have the following two-dimensional 
cross-flow solution 


< p 

T cross 


(x,r,t) 


i iwt 
4 6 


q(x) • Hj 2) 



(5) 


which satisfies the two-dimensional Helmholtz equa- 
tion. 

1 w 2 

0 rr + r^ > r + c 2 " < ^ = 0 (6) 

Obviously, the cross-flow potential equation (5) 
does not represent a complete solution to our problem. 
There is an influence of the body pulsation in front of 
the cross section, x = constant, which has to be taken 
into account. Generalizing a steady -flow concept first 
introduced by K. Oswatitsch [11], we may term this 
part of the flow field "spatial influence" and obtain it 
by subtracting equation (5) from equation (4). 


= -H^ [ r cot a \l(y - a) (y - o'")] ^ > 

o 1 ' ' w n y > a" 

= — K f r cot a *J(y - a) (a" - y)1 r > 01 . 

7T o 1 7 n a" > y 

(7b) 

This expression for the spatial influence gives a 
much clearer insight into the parameters on which the 
flow field depends, namely, on the quantities y, c?, 
and a". 

Since "a and a" represent the wavelengths of the 
advancing and receding acoustic disturbances in the 
x - direction, we see that the solution critically de- 
pends on the ratio of the cylinder deflection y to the 
wave numbers gF and a". In particular, if one of the 
acoustic wavelengths coincides with the deflection 
wavelength, the spatial influence becomes infinite. 
Thus a similar singular behavior of this solution is 
seen to occur as for steady transonic flow. As a 
matter of fact, this breakdown of linearized theory is 
due to the same transonic singularity. This can be 
shown as follows: 

The argument in equation (7a) and (7b) can also 
be written in the form 


Since the integrand in equation (4) decreases 
rapidly with increasing distance of the source-element 
q(£) from the considered station x, we may make an 
approximation by extending the integration from x = - 


showing that the solution depends on the "equivalent" 
steady Mach numbers 
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Since the standing wave pulsation can be obtained 
by superposition of two traveling waves with the veloc- 
ities ± — , our problem resolves into two steady flow 

y 

problems over a wavy walled cylinder with the Mach 

numbers M ± — . The breakdown of linearized theory 
cy 

then occurs when one of these Mach numbers ap- 
proaches one, which corresponds to the previously 
stated coincidence condition. 

IE. STEADY SUPERSONIC FLOW OVER A WAVY 
WALLED CYLINDER 


The foregoing considerations showed that the flow 
over an infinitely long vibrating cylinder can be re- 
duced to two steady flow problems. Therefore, a 
study of the steady wavy walled cylinder flow may give 
considerable insight into the aerodynamic character- 
istics of vibrating cylinders. 

Our objective here is to study the influence of 
wave-amplitude to wavelength ratio to find the effect 
of aerodynamic nonlinearities. This will be done by 
comparing the linearized solutions for an infinitely 
long cylinder with a nonlinear method of characteris- 
tics solution. 

Thus, we assume a sinusoidal wave pattern of the 
form R (x) = R cos yx. 

Specialization of equation (7) to steady flow yields, 
after the proper satisfaction of the boundary condition 
at the body surface, the following expression for the 


pressure coefficient: 

c - - ^ 1 A 1 cos (yx + 9), 

( 8 ) 

P n/m 2 - 1 


where 


I A| = <s/ a 2 + b © = tan -1 — 

a 

(9a) 

J o (X) J,(A) + Y (A) Yj(A) 

a [Ji(X)l*+[T,(X)l I 

(9b) 

J 0 (X) Y,(A) -J,(A) Y q (A) 

(90) 

b " [Ji(A)] 2 + [Y!(A)] 2 

A = y R 's/m 2 - 1' . 
0 

(9d) 


From this linearized solution, two limiting cases 
may be easily obtained by proper expansion* with re- 
spect to the argument A, namely, the slender-body 
solution for values of A « 1 and the Ackeret solution 
for A » 1. By asymptotic expansion for small values 
[ 13] , we obtain the following "slender-body" expres- 
sion for the pressure coefficient: 

C - - 2 y 2 R q R 7(in|+ C ) 2 + j cos (yx + 0 ) 

where 

© = tan " 1 7 . ( 10 ) 


This result could also have been obtained by using 
the well known slender body theory for bodies of revo- 
lution, e. g. , equation 9. 34C Reference 12, and ap- 
plying it to the infinitely long w? yy walled cylinder. 

On the other extreme, by asymptotic expansion 
for large values o f A, we obtain the familiar Ackeret 
solution [ 12, p. 13] 

C = - — sin y x . (11) 

p 

In Figures 3 and 4 the range of validity of these 
two asymptotic solutions is compared with linearized 
theory. Slender-body theory is seen to deviate sharply 
for higher values of A, whereas Ackeret theory turns 
out to be a rather good approximation down to rela- 
tively low values of A. 

In panel flutter calculations exact linearized theory 
or various approximations such as Ackeret theory, 
piston theory, etc. , are used to calculate the general- 
ized aerodynamic forces. Since the flutter amplitudes 
are usually very small, an analysis of the aerodynamic 
forces may well be based upon linearized theory, and 
any deviations of the measured pressure distribution 
should be attributed to boundary layer effects. How- 
ever, some recent experiments by McClure ( 14] and 
Anderson [ 15] , along with currently planned tests of 
wavy wall models in the transonic and low supersonic 
range [ 16] , suggested the desirability of a more ac- 
curate flow field analysis. 

For this purpose a method of characteristics so- 
lution for supersonic flow of an ideal inviscid non-heat- 
conducting gas has been adapted to the present prob- 
lem. This solution, worked out by J. L. Sims and L. 
Saunders* , is based upon the nonlinear non-isentropic 
equation, 


* Aero-Astrodynamics Lab. 


26 



6 


24 - 

A 



FIGURE 3. ABSOLUTE VALUE OF PRESSURE 
COEFFICIENT 



FIGURE 5. WAVY WALL CYLINDER INVESTI- 
GATED BY ANDERSON IN HIS WIND-TUNNEL 
TESTS 
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FIGURE 4. PHASE ANGLE OF PRESSURE 
COEFFICIENT 


Phase Angle 0 (degrees) 

2 n„ ,,2 n«r 2uv 9u 


u“ 8u . v‘ 9v 


9r 




uv 9s 
vVR 9n * 


( 12 ) 


where x, r, and 0 are cylindrical coordinates, u, v, 
w are the velocity components of the velocity V, c the 
speed of sound, R the gas constant, s the entropy and 
n the direction normal to the streamline. 


We will apply this numerical method of charac- 
teristics solution to a wavy walled cylinder whose di- 
mensions are shown in Figure 5. This model has been 
extensively studied by Anderson [ 15] , and we refer to 
his report for details of his experimental investiga- 
tions. Our purpose here is to show a comparison 


FIGURE 5a. SKETCH SHOWING DIMENSIONS OF 
ONE WAVELENGTH OF THE ANDERSON MODEL 


between nonlinear theory for a finite length model and 
linear theory for an infinitely long cylinder. Figure 
6 shows the variation of the pressure coefficient for 
supersonic flow at M=3. 0 over the first seven waves. 
Linear theory, of course, is a pure sine wave. Non- 
linear theory also exhibits the same sinusoidal pattern 
and, in particular, shows no leading edge effect. The 
only effect here of aerodynamic nonlinearity is an 
overall shift in pressure levels which causes an inte- 
grated compressive loading of the shell. A similar 
comparison for M= 4. 62 shows the same general char- 
acteristics (Figure 7). However, nonlinear theory 
begins to show a slight leading edge effect and a no- 
ticeable deviation from the sinusoidal variation. Fur- 
thermore, a small upstream phase shift of the peak 
amplitude is to be observed. These characteristics 
[are further enhanced for flow at M = 10. 0 as shown in 
Figure 8. Here a quite significant leading edge effect 
is seen to occur. 


In exploring the lower supersonic Mach number 
range, computer runs for M = 2. 0 and M = 1. 5 show 
essentially the same results as for M = 3. 0. A further 
extension to the range 1 < M < 's/I!? however, could not 
be obtained for the present configuration (Figure 5) 
because the flow is soon decelerated to local subsonic 
Mach numbers in going through the first shock wave. 
This prohibits further analysis by the method of char- 
acteristics. Only a sharp decrease in wave amplitude 
allow ed a further analysis in this Mach number range. 
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FIGURE 6. PRESSURE DISTRIBUTION ALONG WAVY WALLED CYLINDER OF FIGURE 5 

M = 3. 0 a/L = 0. 025 



FIGURE 7. PRESSURE DISTRIBUTION ALONG WAVY WALLED CYLINDER OF FIGURE 5 

M = 4. 62 a/L => 0. 025 
00 
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FIGURE 8. PRESSURE DISTRIBUTION ALONG WAVY WALLED CYLINDER OF FIGURE 5 

NLo~ 10, 0 a/L -- 0.025 


IV. STEADY SUPERSONIC FLOW OVER A TWO- 
DIMENSIONAL WAVY WALL 


The foregoing results for the wavy cylinder flow 
shall now be further substantiated by an analysis of the 
two-dimensional wavy wail flow. Here again the same 
basic features are found as for the cylinder flow, 
namely, no leading edge effect for medium Mach num- 
bers and an overall shift in pressure levels. As an 
example, Figure 9 shows the pressure distribution 
along a two-dimensional wavy wall at M - 1.6 and an 
amplitude to wave length ratio a/L = 0. 02 by linear 
and nonlinear theory. For this Mach number, the non- 
linear method of characteristics solution predicts 
supersonic flow throughout the whole flow field. Re- 
ducing the free-stream Mach numbei to M = 1. 405, 
however, already leads to subsonic flow right after the 


seventh shock (Figux^e 10). Lines oi constant Mach 
number for this particular case are shown in Figure 
11. A similar analysis for an amplitude to wave length 
ratio a/L = 0.01, on the other hand, reestablishes 
supersonic flow throughout the whole flow field at 
M- 1.405. To illustrate this amplitude dependence 
of the flow field a little further, Figures 12 and 13 
show a thiid example, anamplitude to wavelength ratio 
a/L = 0. 04. At this higher amplitude, the flow be- 
comes subsonic already at M - 1. 6 in going through 
the second shock wave. 

McClure [ 14] gives measurements of the pressure 
distribution over two-dimensional wavy walls with am- 
plitude to wavelength ratios of 0. 02 and 0. 01 at 
M = 1.405. The present calculations may point out 
the appreciable effects of aerodynamic nonlinearities 
in the low supersonic flow regime. The available ex- 
perimental and analytical information in this particular 



FIGURE 9. PRESSURE DISTRIBUTION ALONG A TWO-DIMENSIONAL WAVY WALL 

M =1.60 a/L = 0.02 

00 
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speed range is unfortunately not sufficient to evaluate 
this effect in more detail. It te, therefore, one ob- 
jective of the currently planned experimental program 
[ 16] to investigate these flow aspects ver" carefully. 


V. CONCLUSIONS 


a/I. -.04 M -1.6 

«c 



FIGURE 12. NONLINEAR FLOW FIELD OVER 
TWO-DIMENSIONAL WAVY WALL 
a/L = 0.04 M = 1. 6 

CO 


Recent panel flutter studies [7, 8, 14] have led to 
the conclusion that linearized potential flow theory is 
inadequate for an accurate description of the flutter 
phenomenon at the low supersonic Mach numbers. 
Two phenomena were suggested as the main reasons 
for the apparent disagreement between theory and ex- 
periment, i. e. , the effect of aerodynamic nonlinearity 
and the boundary layer effect. S<nce the flutter dis- 
placements are usually very small, it seems plausible 
to ignore the nonlinear transonic flow effects and to 
ascribe the dominant effect to the neglect of viscosity. 
Recent studies by McClure [ 14] and Anderson [ 15, 17] 
attempted to predict this boundary layer influence and 
to verify it by measuring the flow over wavy surfaces. 

For the purpose of this study, the steady super- 
sonic flow over wavy walled planai and cylindrical 
configurations was analyzed by means of a nonlinear 
method of characteristics solution and compr-rad with 
linearized solutions. There is an appreciable depen- 
dence of the flow characteristics upon the amplitude to 
wavelength ratio causing mixed flow phenomena for 
low supersonic free-staeam Mach numbers. 


a/L - . 04 M *1.6 




A more detailed anal^-ds using transonic flow 
theories is necessary in order to evaluate the signifi- 
cance of transonic nonlinearities upon panel flutter. 
Conceivably, large amplitude flutter, flutter of buckled 
panels, and flutter of slightly curved panels could be 
appreciably influenced by such transonic flow effects. 
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III. AEROTHERMODYNAMICS 



CONTINUOUS CARBON ABSORPTION STUDIES 
By 

C, C. Ferriso* , C. B. Ludwig* , and Robert M. Huffaker 




SUMMARY 


For thick optical paths of interest for Saturn type 
1st stage boosters, carbon contributes more than 50 
percent of the total emission. The absorption co- 
efficient of this carbon produced in a RP-i/oxygen 
rocket engine combustion chamber has been experi- 
mentally determined and an analytical model developed 
which accurately expresses the absorption coefficient 
as a function of temperature, density, and wavelength. 
These measurements were made on a small rocket 
motor with the perfect nozzle designed by the method 
of Foelsch, This nozzle produces homogeneous, 
shock -free, axially parallel exit gases. 


I. INTRODUCTION 


The radiative heating of the base of large booster 
vehicles using a hydrocarbon-oxygen propellant is pri- 
marily due to emission from three constituents of the 
exhaust gases*; C0 2 , H 2 0, and solid carbon. The emis- 
sion from the molecules is restricted to the various 
vibration rotation bands and has a complex dependency 
on frequency. The solid carbon, on the other hand, 
gives rise to a continuous and smoothly varying spec- 
trum. Water and carbon dioxide are equilibrium prod- 
ucts of the combustion; solid carbon is not. Thermo- 
dynamic calculations predict no solid carbon for the 
mixture ratios at which the engine is normally oper- 
ated. 

An estimate of the relative importance of the dif- 
ferent radiating sneeies may be obtained from an 
examination of available data for the total hemispheri- 
cal emissivities for conditions characteristic of the 
exhaust. It is found that, for the thick optical paths of 
interest for Saturn type boosters , carbon contributes 
more than 50 percent of the total emission, The data 
for solid carbon used in this estimate were obtained 
from theoretical values for the absorption coefficients 
(Fig. i) ofcarbon particle clouds which, in turn, were 
based on the properties of pure bulk graphite. The 
experimental determination of these absorption coef- 
ficients for the carbon produced in RP-i/oxygen rocket 
engine combustion chambers is the primary purpose 
of this work. 

* Employed by General Dynamic^/Convair 
San Diego, California 


The determination of the emissive properties of 
volumes of gas containing dispersed carbon particles 
is, in some aspects, much simpler than the equivalent 
problem for the molecular emitters, and in some as- 
pects much more difficult. The molecular absorption 
spectra, when observed at high spectral resolution, 
generally show a highly complex functional dependence 
on wavelength. In contrast to this the carbon absorp- 
tion spectra are expected to vary smoothly with wave- 
length. Thus Beer's law applies to- transmission or 
emission measurements made at moderately low spec- 
tral resolution and measurements made on thin sam- 
ples of gas may therefore be extrapolated reliably to 
any desired thickness. This greatly simplifies the ex- 
perimental problem since only the temperature and not 
the pathlength need be varied. 

Since carbon is present in particulate form, the 
particles will scatter, as well as absorb light. If the 
particles are not small compared to wavelength, this 
scattering may significantly affect the radiant heat 
transfer. In general, the effect of scattering will be 
to lower the total energy radiated by the gas from that 
expected when scattering is neglected, although the 
radiation in particular directions may actually be en- 
hanced at the expense of that in other directions for 
some fiow field configurations. 

One of the most difficult aspects of the present 
problem is the establishment of the relation between 
the amount and properties of the carbon particles pro- 
duced in the small scale rocket engines with those of 
the full scale engines. The solid carbon produced in 
the combustion of a fuel -rich mixture and oxygen is a 
nonequilibrium product. In the absence of reliable 
quantitative methods for predicting the amount formed 
in a given situation, some measurements on the full 
scale engine will be required. 


IL SCOPE OF THE PRESENT EXPERIMENTAL 
PROGRAM 


The experiments consist of measuring the emis- 
sion from small, hydrocarbon -fueled rocket engines 
operated in combustion regimes where carbon is 
formed. These measurements are made at the exit 
plane using the techniques and equipment applied to 
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FIGURE 1. LINEAR ABSORPTION COEFFICIENT 
FOR CARBON PARTICLE CLOUDS 


measurements of emissivities of gases. By holding 
combustion constant and varying the expansion ratio 
of the engine, one can vary the temperature of the 
system while holding ‘he mass fraction of carbon 
nearly constant. The d . & are reduced to the form of 
absorption coefficient as a function of combustion con- 
ditions: temperature and wavelength. Since the carbon 
is produced by an essentially unknown process in the 
combustion chamber, the application of data taken in 
small scale experiments to the full scale motor must 
be carried out with some care. Although it is not pos- 
sible to duplicate the combustion conditions of the large 
engine in the small scale engine, part of the program 
is devoted to examining the effect on carbon content 
and the spectral distribution of the absorption coeffi- 
cient of varying different combustion parameters . The 
object is to identify those design parameters that con- 
trol the carbon formation, such as mixture ratio, 
chamber pressure, chamber length, injector pattern, 
injector hold size, throat diameter, state of oxidizer 
(LOX or GOX) , and cooling water flow rate. 


m. TEST APPARATUS AND MEASUREMENT 


A. ROCKET BURNER 

Absorption coefficient measurements require 
homogeneous samples of gas of known temperature, 
pressure, and composition, containing measurable 
quantities of carbon particles. These gas -particle 
systems are produced by small rocket motors with 
contoured nozzles burning RP-1 and oxygen as pro- 
pellants. 


An exploded view of a fypical motor is shown in 
Figure 2. The motor is rated at 160 lb (715 N) nom- 
inal thrust (0.6 lb/ sec [0.27 kg/sec]mass flow rate). 
The motor consists of an electrolytic copper injector, 
a machined electrolytic copper combustion chamber - 
nozzle piece (the "liner") , and a stainless steel cool- 
ing jacket. 



FIGURE 2. EXPLODED VIEW OF MOTOR 

The liner can be made with different area ratios. 
All liners are constructed with the same chamber di- 
mensions and nozzle contour to the throat. The diverg- 
ing sections are perfect nozzles designed by the meth- 
od of Foelsch, using an IBM 7090 computer program; 
this computes the nozzle contour for a given area ra- 
tio, throat diameter, and combustion product ratio of 
specific heats. These nozzles produce homogeneous, 
shock-free, axially parallel exit gases, as shadow- 
graph studies have shown. Figure 3 shows some of 
the shadowgraph results. 

B. INFRARED EMISSIVITY APPARATUS 

Spectral emissivify measurements are made 
of the carbon continuum between 0. 7 ju and 4. 0 n using 
spectrometer equipment. (Fig. 4) . 

The energy emitted by the burner exit gases is 
focused on the entrance slit of Perkin-Elmer Model 
99 double-pass monochromator. A thermocouple de- 
detector or a cooled PbS detector , and a CaF 2 prism 
is used in the monochromator. A wave number cali- 
bration of the monochromator is made using known 
absorption bands of various gases. For the energy 
calibration, the two-position diagonal mirror is ro- 
tated to fill the monochromator aperture with radia- 
tion from a Barnes Engineering Model RS6A 1000° C 
( 1270°K) blackbody calibration source, A tempera- 
ture measurement is made of the source during cali- 
bration with a Pt-Rh thermocouple placed on the cavity 
surface. The optical system is enclosed and flashed 
with dry nitrogen to reduce atmospheric absorption. 
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FIGURE 3. SHADOWGRAPHS OF BURNER EXITS 



FIGURE 4. SPECTROMETER FOR SPECTRAL 
EMISSIVITY MEASUREMENTS 

C. TEMPERATURE MEASUREMENTS 

Temperatures, both in the continuum and 
molecular spectral regions, are measured with an 


absorption-emission technique. The spectral emissi- 
viiy, e(u), is calculated from the measured trans- 
mission with 


e(p) =1 - 


JM. 


(i) 


where I ( v ) is the incident globar energy and I(^)is 
the transmitted energy of the source (globar) attenu- 
ated by the hot exhaust gas. The temperature can then 
be determined from the corresponding Planck black- 
body function. 


D. SPECTRAL EMISSIVITIES MEASUREMENTS 


The spectral emissivities are determined 
from emission measurements and the experimental gas 
temperature measurements. The energy emitted by 
the exhaust gas is obtained on a calibrated energy ba- 
sis by relating the recorder deflection to the blackbody 
calibration source. Readings from the recorder traces 
(exhuustgas and calibration source) are made at con- 
stant monochromator drum settings and serve as input 
to a computer data reduction program that converts 
drum settings into wave numbers by using the original 
wave number calibration curve. The computer pro- 
gram then evaluates the spectral emissivities at a par- 
ticular wave number, r, from 


G(v)N*(v # T b ) 
B~(v) N°( ^T) 


(2) 


where G(v ) and B(y) are the voltage deflections at a 
particular wave number for the hot gas and blackbody 
source, respectively; N°(^, T 3 ) is the Planck black- 
body energy value for the calibration source at its 
temperature, T^; and N°(y,T) is the Planck black- 
body energy value for the hot gas where T at the meas- 
ured mixture ratio is taken from a curve fit to the ex- 
perimental temperature data points. The data reduc- 
tion program is also set up to give total emissivities, 
spectral absorption coefficients, and simple radiance 
(energy) spectra. 


IV. RESULTS 


Data have been collected on two liners having area 
ratios of 5. 25:1 and 3:1. The exit temperatures ver- 
sus mixture ratios are plotted in Figure 5, together 
with theoretical values, assuming complete equilib- 
rium in the chamber and shifting equilibrium expansion 
to the exit plane, which represents an upper limit. It 
is interesting to note that the experimental values ex- 
ceed the theoretical predicted values. It is believed 
that this behavior is, in part, due to the non -equilibrium 
process of condensing carbon in the chamber. 
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FIGURE 5. PLOT OF EXPERIMENTAL AND 
THEORY VALUES OF EXIT TEMPERATURES VS 
MIXTURE RATIOS FOR THE 3:1 AND 5.25:1 
ROCKET MOTORS 

In Figures 6 and 7 we show typical traces of the 
measured spectral emissivities at a mixture ratio of 
1.42 for the 5. 25: land the 3:1 area ratio number, re- 
spectively. The molecular emission at 4. 3 p and 2.7 p 
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FIGURE 6. TYPICAL TRACE OF THE MEASURED 
SPECTRAL EMISSIVITIES AT A MIXTURE RATIO 
OF 1.42 FOR THE 5.25:1 AREA RATIO MOTOR 

is clearly indicated. These traces are reduced to pure 
continuum radiation by subtracting the molecular radi- 
ation from C0 2 and H?0. In Figure 8 is plotted the 
"effective spectral absorption coefficient" of carbon 
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FIGURE 7. TYPICAL TRACE OF THE MEASURED 
SPECTRAL EMISSIVITIES AT A MIXTURE RATIO 
OF 1.42 FOR THE 3:1 AREA RATIO MOTOR 
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FIGURE 8. COMPARISON OF EXPERIMENTAL 
AND THEORETICAL DATA WHERE THEORETICAL 
CALCULATIONS ARE ADJUSTED AT 2.2 MICRONS 

times the carbon density, p c » as dots. The spread 
indicates the experimental uncertainty, as well as the 
uncertainty in the data reduction procedures. The 
dashed lines represent the model calculations taken 
from Figure i. These lines were placed to give the 
best agreement with the experimental points. There 
is good agreement in the shape of these curves. The 
scaling factor is the unknown carbon concentration, 
using the theoretical absorption coefficients. The 
carbon concentration is given in Figure 9 for all mix- 
ture ratios investigated. The agreement between the 
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«• 

values of the two different engines indicates the proper It is found that the carbon concentration is strongly 
temperature dependency of the absorption coefficients. dependent on the mixture ratio. 



MIXTURE RATIO 


FIGURE 9. CARLDN CONCENTRATION FOR ALL 
MIXTURE RATIOS INVESTIGATED 
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INFRARED SPECTRAL TRAMSMITTANCES FOR DETERMINING RADIANT FLUX 

OF INHOMOGENEOUS GASES 

By 

Robert M. Huffaker 




SUMMARY 

^ -r ? ^ -T" 


A general method has been developed to calculate 
spectral transmittance s of inhomogeneous gases from 
the properties of homogeneous gases. Thus, any cal- 
culation of jpectral transmitt?.nce fo v a particular in- 
homogeneous gas path can be done by properly com- 
bining known "band model" data on gases at constant 
temperature, pressure, and concentration. The meth- 
od developed is based principally on two special spec- 
troscopic concepts: (1) the molecular band -model, 
and (2) the Curtis-Godson approximation. This 
method has been experimentally ver ' p ied over the tem- 
perature, pressure, and wavelength range typical of 
Saturn V exhausts. ^ 


I. INTRODUCTION 


The exhaust plume from a multiple -engine con- 
figuration like the Saturn-V is an extraordinarily 
complex gas-dynamical structure, „ The exhaust jets 
from the different engines interact with the atmosphere 
and with each other through an intricate series of shock 
waves and turbulent mixing layers. Some of the ex- 
haust material flows back into the base region. Where - 
ever the exhaust mater ial mixes ith the atmosphere, 
it can burn, and, at low altitudes, it dees. With such a 
complicated flow field controlling the heat tranjfer to 
a base region, it is not surprising that simple models 
of the heat transfer must be tied strongly to empirical 
data. 

In an effort to learn more about the base heating 
problem and to provide both mathematical techniques 
and values of parameters required for solution of the 
problem, a number of studies have been undertaken, 
These involve descriptions of the flow field, the radi- 
ative properties of its constituents, and a radiation 
calculation method for inhomogeneous gases. 


Heat transfer to the base of a launch vehicle can 
occur through two modes: convection and radiation. 
There is evidence that, on large vehicles burning RP-i 
(kerosene) and oxygen, these two modes are both im- 
portant. There currently exists no satisfactory tech- 
nique of calculating base heating for a given engine 
array at various altitudes and vehicle velocities. We 
therefore require both better means of measuring heat 
transfer to the base region of vehicles actually under 
test, and more rigorous and reliable means of pre- 
dicting heat transfer to vehicles before testing. 

Convection to the base region depends upon the 
properties of the gas flow close to the surface being 
heated. Simple models c' ** • low, together with 
boundary layer correlations lor the heat transfer co- 
efficient, have been used in the past to arrive at local 
heating rates on the heat shield. Radiation, on the 
other hand, is sensitive to the entire flow field not 
occluded by the vehicle structure. Here predictions 
have been made either by representing the exhaust 
plume as a surface with a given emissivity, tempera- 
ture, and configuration, or by assuming a simple geo- 
metrical configuration for the het gas, which is pre- 
sumed to be uniform , and integrating the equation of 
transfer with a simplified absorption coefficient. 


Water vapor, carbon dioxide, anu solid carbon 
particles (soot) are the principal radiators in the 
firsi -stage plume. Water vapor is the only radiator 
of importance in the plumes of the hyde*. gen-fueled 
stages. Extensive experimental and theoretical stud- 
ies of the emission of COg under various conditions 
have been made, and the properties of this molecule 
are v/ell understood. Much less is presently known 
about H 2 O. A NASA sponsored program, directed by 
Dr, Carmine Ferriso, is now underway at General 
Dynamies/Convair. The program was set up to meas- 
ure the emission and absorption chax'acteristics of 
water vapor for temperatures from 300°K to 3000°K 
and for pressures of 0.1 atmospheres to 1.0 atmos- 
pheres. Band model parameters will also be obtained 
for water vapor from this experimental data and addi- 
tional theoretical studies. 

To calculate the radiant heating to a rocket base 
by the exhaust plumes, one must know the spectral 
transmittance s of the plume gases. Since plumes are 
strongly inhomogeneous, one must determine spectral 
transmittance s of inhomogeneous gases explicitly to 
develop a calculi >n method. A general method has 
now been developed with the Warner and Swasey Com- 
pany Contract under the direction of Dr. Richard 
Tourin, to calculate spectral transmittances of 
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inhomogeneous gases from the properties of homo- 
geneous gases. Thus, any calculation of spectral 
transmittance for a particular inhomogeneous gas path 
can be done by properly combining known data on gases 
at constant temperature, pressure, and concentration. 

The method developed is based principally on two 
special spectroscopic concepts: (1) the molecular 
band model, and (2) the Curtis -Godson approximation. 
The molecular band model is a mathematical repre- 
sentation of the effective radiant energy in portions of 
the spectrum where gas molecules emit and absorb 
energy. The band model yields an explicit, closed 
formula for the molecular radiation within each se- 
lected spectral region of interest, which uses as input 
data the averaged line strength, the averaged line spac- 
ing and the averaged line half width. An average of 25 
wave numbers is considered sufficient. An exponen- 
tial line intensity distribution or a delta function dis- 
tribution is also used in the band model. The use of 
a band model is critically important for practical cal- 
culations of gas radiation. Without a band model, one 
would have to deal individually with each of the many 
thousands of spectral lines that contribute to molecular 
radiation, a much bigger and more difficult problem, 
parts of which are beyond solution at present. The 
Curtis -Godson approximation is a method of combining 
the parameters that appear in the band model formulae 
in such a way that the parameters needed for an in- 
homogeneous gas calculation are obtained solely from 
homogeneous gas data. This is also a critical factor, 
because the necessary parameters need be calculated 
or measured in the laboratory only for uniform gas 
samples. Without the Curtis -Godson, or some equal- 
ly good approximation, one would have to treat each 
inhomogeneous gas path as a special case. 


II. THEORY 


mittanceof the inhomogeneous gas. The method used 
to calculate radiant flux for a homogeneous gas will be 
outlined first; this should clear the approach to the in- 
homogeneous case. 



FIGURE i. SCHEMATIC OF ROCKET EXHAUST 
PLUMES. THE ARROW REPRESENTS A TYPICAL 
INHOMOGENEOUS GAS PATH 

Assume a mass of homogeneous hot gas, which 
emits and absorbs radiant energy at wave number v , 
A beam of radiant energy afc wave number v 9 of arbi- 
trary magnitude I q (v) is sent through the hot gas. 
Some of this beam is absorbed in the gas, and the re- 
mainder I(v) emerges. The gas itself emits radiant 
energy at wave number v % of magnitude N (v,T). 

O 


A. PROBLEM STATEMENT AND DEFINITION 

We wish to predict the radiant energy trans- 
ferred from an inhomogeneous mass of hot gases. The 
radiated energy is distributed irregularly over a broad 
range of the spectrum, mainly at wave numbers (the 
frequency divided by the velocity of light) between 2000 
cm -1 and 10,000 cm -1 . Once the energy is obtained 
along a line of sight, the radiant flux is then integrated 
over all spectral, regions for this line of sight, and 
integrated over solid angle to obtain the total radiant 
flux. In this report we are concerned with the first 
part: how to determine the radiant flux incident upon 
the base in any spectral range , contributed by the hot 
gases along a line of sight like the one in Figure i. 
The key variable turns out to be the spectral trans- 


The spectral radiance of the hot gas Ng( v, T) is 
defined as the radiant power per unit area, per unit 
solid angle, per unit wave number interval, in watt- 
cm’ 2 -steradian -1 -cm, . 


The spectral transmittance t ^ ^ is the fraction of 
incident radiant energy emerging from the gas; i.e.. , 


t(v) = 


Uv) 

V”) 


( 1 ) 


The spectral absorptance a(v) is the fraction of 
incident radiant energy absorbed in the gas; i,e. , the 
fraction not transmitted, or 


1 -t(v) = 1 -I(m)/I o (i>). (2) 
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The spectral emissivily e(y) is the ratio of the 
gds radiance, N (y,T) to the radiance of a blackbody 
Kirchhoff's law states that 

e(v) = a(v) , (3) 

from which one can easily get e(v) by measuring a(v). 
In radiant transfer studies of isothermal gases, atten- 
tion is focused on the spectral emissiviiy , because if 
one knows e(y) and the temperature, the radiance is 
easily computed; thus, 


N (g.T) = c(v,T) N fa (u,T) . (•-) 

hi considering strongly inhomogeneous gases, the 
radiance calculation is formulated quite differently, 
because the emission and absorption within the gas 
must be considered explicitly. The equation of radi- 
ative transfer becomes 


n 


model is a mathematical formula which gives the 
transmittances, for intervals Ay in which one is in- 
terested, without the need to deal with the individual 
spectral lines. The band model is vitally necessary 



FIGURE 2. SCHEMATIC INFRARED SPECTRUM 
OF A HOT GAS IN TWO NARROW INTERVALS, 
PLOTTED AS TRANSMITTANCE VERSUS 
WAVE NUMBER 


where the terms are understood by reference to Figure 
1. One now considers the inhomogeneous structure of 
the line of sight, and considers it subdivided into a 
series of zones, each of which is homogeneous. The 
quantities in equation (5) are then; H(y i9 Ay.) is the 

irradiance at the point of interest within wave number 
Ayj centered at wave number yj; N^VpT.) is 
blackbody spectral radiance at yj, at the temperature 
Tj of the ith zone; r. j(y^Ay^) is the spectral trans- 
mittance of the inhomogeneous path comprising all 
zones between the ith zone and the point of interest; 

and T,(y.,Ay.) is the spectral transmittance of the 
J 1 

path including the ith zone and all other zones between 
the ith zone and the target point. The Planck radiance 
functions are well known, and tabulated, so that the 
quantities we seek are the inhomogeneous transmit- 
tances; Warner and Swasey Company has accomplished 
this by applying the Curtis-Godson approximation. 
First, the band model method for a homogeneous gas 
will be defined. Thenitwillbe shown how these homo- 
geneous gas parameters are combined to get the cor- 
responding parameters for the inhomogeneous case. 


B. THE STATISTICAL BAND MODEL 


in radiant heat transfer calculations for two reasons: 
(i) the number of spectral lines involved in a prac- 
tical problem may run into hundreds of thousands, and 
for most of these not enough is known for accurate 
calculations; and (2) practical spectroscopic instru- 
ments distort the spectral lines, so that the "true" 
spectrum cannot be measured. The statistical band 
model formula for the spectral transmittance r at wave 
number v is 

In £ = 2 n ■£ f(x) , (6) 

where 



The term is the ratio of average line width to aver - 

age line spacing in the spectral interval Ay at v; — is 

the ratio of average line strength to average line spac- 
ingin the same interval; f(x) is the Ladenberg-Reiche 
function; and i is the path length in the gas. The 

g 

strength parameter t is proportional to P , the par- 
cl a 

tial pressure of absorbing gas; thus, 


While many band models have been studied, 
the statistical model best represents the spectra of 
water vapor and carbon dioxide. The need for band 
models is illustrated in Figure 2 where two small spec- 
tral intervals Aiq and Ay 2 are illustrated showing typ- 
ical spectral lines of a molecular spectrum, A band 


S _ S° 
d " d a * 

Then, if we define 

M = P I, 
a 


41 


equation (7) can be rewritten as 


x 


S- 

d 


M 


2 7T 


T 

d 


(7a) 


To apply the statistical model to calculating trans- 
mittances of homogeneous gases, we need only two 
parameters: the strength parameter S/d, and the 
width parameter y/d, averaged over each spectral in- 
terval of interest. This information is being obtained 
under contract to General Dynamics/Convair. 

C. COMBINATION OF ISOTHERMAL BAND 
MODEL PARAMETERS TO GET INHOMO- 
GENEOUS TRANSMITTANCES 

Several techniques were investigated by 
Warner and Swasey for combining isothermal band 
model parameters to predict the transmittance through 
an inhomogeneous path. The Curtis -Godson approxi- 
mation gave the best results. 

The Cutris-Godson method is based on the pre- 
mise that one can substitute a hypothetical homo- 
geneous gas path that would have the same spectral 
path will have the same spectral transmittance as an 
inhomogeneous path if the parameters (y/d)' andX' 
as defined by equations (8; and (9) below are used in 
the band model formulas for the transmittance of a 
homogeneous gas, equations (6) - (7a). 




l 

h-t 


h . 8 \ 
P„ < 7 ) 


i/a\ h X hT 

L<r> v ( r> p bJ 


V ( ff 

ft d 


( 8 ) 




) P. 
h 


h + S 



( 9 ) 


The subscript h in equations (8) and (9) stands for 
each of the homogeneous zones into which the inhomo- 
geneous path has been divided. The limits indicate 
that the number of zones in the inhomogeneous path is 
i. Substituting equations (8) and (9) into equatipns 
(6) and (7) will give Tj for use in equation (5). When 
the limit 1-i is used instead of i, one can calculate 
of equation (6). P, is the pressure of nonradiat- 

Y b 

ing gas in each zone, such as nitrogen, and (~p) is 


the line -width parameter for this gas. The derivation 
of equations (8) and (9) is given in G. R. Harrison, 
Practical Spectroscopy, Chapter 14. Equations (8) 
and (9) are derived by requiring that the mean trans- 
mittance of the inhomogeneous path is exactly equal to 
the mean transmittance of the homogeneous path for 
the extreme cases of infinitesimally small path lengths 
and very large path lengths. For intermediate path 
lengths, the equality is only approximate. 

D. RECAPITULATION 

The Curtis-Godson method is based on the pre- 
mise that one can substitute a hypothetical homogene- 
ous gas path that would have the same spectral trans- 
mittance as a real, inhomogeneous path encountered 
in a rocket exhaust plume. These parameters were 
shown to be calculable from the band model parame- 
ters of the homogeneous zones. Once the parameters 
of the hypothetical path were determined, the spectral 
transmittance was calculable from the statistical band 
model formula for spectral transmittance. The com- 
plex problem of determining transmittances of inho- 
mogeneous gases has thus been reduced to the simpler 
problem of determining band model parameters of ho- 
mogeneous gases. 


HI. EXPERIMENTS AND RESULTS 


Experiments were carried out at Warner and 
Swasey to test the theory, determine the accuracy, 
and obtain values of band model parameters for spe- 
cific cases. Band model parameters were measured 
for isothermal gases at various temperatures, the re- 
sults were used to predict transmittances of known 
inhomogeneous paths, and the inhomogeneous path 
transmittances were compared to the theoretical pre- 
diction. 

Using a furnace -gas cell arrangement, numerous 
transmittances of water vapor inhomogeneous paths 
were measured and compared to transmittances cal- 
culated by the Curtis-Godson combination method in 
conjunction with the statistical band model. Each in- 
homogeneous path consisted of two zones, of different 
temperature, pressure, and composition. Zone path 
lengths varied from 1 1/2 inches to 24 inches, and 
total paths from 3 inches to 48 inches. Typical results 
are given in Table 1. Table 1 also shows two-zone 
transmittances calculated from the transmittances of 
the individual zones by the commonly used Beer- 
Lambert law, which(gave’.an order of magnitude larger 
error than did the Curtis-Godson method. 
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FIGURE 3. SCHEMATIC DIAGRAM OF FURNACE- 
GAS CELL-SPECTROMETER SYSTEM INCORPO- 
RATING TWO FURNACES. OPTICAL BEAM 
PASSES THROUGH EACH FURNACE THREE 
TIMES 
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TABLE I 


Two Zone .ter Vapor Transmittance s 
v = 3990 cm 1 


Tj = 637 °K 
i i = 24 inches 


T 2 = i273°K. 
i 2 =24 inches 



Zone 1 

Zone 2 


|two Zone Transmittances i 






Measured 


Calculated 

Case 

H 2 0 

Zonal 

H 2 0 

Zonal 




No. 

Pressure 

Trans. 

Pressure 

Trans. 


Curtis- 

Beer- 


(mm Hg) 

Tl 

(mm Hg) 

T'i 

t 12 

Godson 

Method 

Lambert 

Law 

1 

53 

.711 

49 

.747 


.631 

.531 

2 

53 

.711 

49 

.745 


.630 

.530 

3 

51 

.722 

101 

.573 

.501 

.518 

.414 

4 

51 

.714 

105 

.563 

.492 

.507 

.402 

5 

50 ’ 

.720 

148 

.460 

,413 

.424 

.331 

6 

52 

.710 

148 

.452 

.402 

.415 

.321 

7 

104 

.503 

100 

.597 

.398 

.414 

.300 

8 

100 

.519 

106 

.564 

.384 

.409 

.293 

9 

112 

.474 

102 

.570 

.363 

.383 

.270 

10 

102 

.508 

147 

.459 

.327 

.348 

.233 

11 

101 

.520 

154 

.436 

.317 

.339 

.227 

12 

138 

.395 

138 

.474 

.276 

.294 

.187 

13 

151 

.358 

151 

.443 

.247 

.260 

.159 

14 

147 

.373 

153 

.437 

.243 

,267 

.163 

15 

152 [ 

.358 

154 

.436 

.240 

.258 

.156 
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IV. CONTROL STUDIES 



SOME SUFFICIENT CONDITIONS FOR MINIMAX CONTROL 



By 

Thomas E. Carter 


SUMMARY 




^ Presented are the initial efforts of a mathematical 
investigation to determine sufficient conditions for so- 
lution of a class of optimal control problems in which 
the quantity to be minimized is the peak value of some 
measure of performance of the control. This work 
was motivated by the fact that many guidance and con- 
trol problems can be formulated as problems of this 
class, yet very little theoretical work has been done 
in this area. Some theorems are developed which ap- 
ply to this class of problems. A simplified version of 
a booster control problem known as Ihe !, constrained 
bending moment minimax drift problem" is used as a 
model for application of the theorems. Some conjec- 
tures are made regarding extending the work to other 
minimax problems. 




a , 


I. INTRODUCTION 


A class of problems which has received very little 
attention in the literature of optimal control theory is 
that in which the function to be minimized is given by 
the peak value of a function over a time interval. 
Problems such as this belong to a class which liave 
been called "minimax" problems atMSFC. These in- 
clude 'guidance and control problems of space vehicles 
in which it is desirable to minimize peak values of some 
selected index of vehicle performance such as lateral 
drift or bending moment. 


booster control problem which has been under consid- 
eration at MSFC for some time. This problem is 
known as the constrained bending moment minimax 
drift problem, A minimax control law is developed 
for a simplified version of this problem. 


II. DEFINITIONS AND STATEMENT OF THE 
PROBLEMS 


Consider the following set of equations , 

J? = f(x, t, u (t) ) , (1) 

defined over a closed interval I C R where R is the 
set of real numbers. Differentiation is with respect 
to t e Ip. The real valued function u(t) belongs to a 
class U called the set of "admissible controls." The 
state vector x belongs to region X called an "admis- 
sible set" of n-vectors. For eachx e X and u(t) e 
U and subinterval IC I F , it is assumed that only one 
solution x(* o ,t;u(t)) exists such that x(x q , t Q ; u(t)) 

= x . 
o 

Associated with (1) is a scalar y (t,u( t) ) which 
denotes the "cost" of u(t) for each x , This is given 
by the equation 

y = h(t, x(x Q ,t; u(t)), u(t)). (2) 

For a given x q e X this may be written as 


Since control theory so far has had so little to say 
about cost functions of this type, it was necessary to 
go back to the basic techniques of real analysis to in- 
vestigate this type of problem. By use of this ap- 
proach, sever l theorems have been developed which 
are applicable to minimax type problems. These 
theorems can be used to find a solution to some simple 
minimax problems, but no systematic method can yet 
be claimed. It is hoped that other theorems can be 
added to those in this report so that a general theory 
of minimax control elements can be built up similar to 
that which exists in other areas of optimal control. 


y = g(t, u (t) ) , (3) 

For y(t Q , u(t) ) - y , we have 
t 

y(t, u(t) ) = f g(T, u(r))dr+ y Q . (4) 

*o 

Minimax elements with respect to (1) and (3) are 
defined as follows: 


This report gives a mathematical statement of the 
problems including definitions of minimax elements 
and supreme minimax elements. The general theo- 
rems are then developed and applied to a simplified 


Consider a closed interval IC I and a point 
t Q € I such that X(x q , t o ;u(t)) = x Q and y(t Q ; u( t) ) = y . 

If there exist a control u'(t, I) e U and a point y'e R 
with the property that 
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y' = min max.y(t,u(t)) , 
ueU tel 

then 

(1) y'(I) is called a minimax value ofy(t, u(t)) 
with respect to I, 

(2) u'(t,I) is called a minimax control, and 

(3) x (x , t; u’(t,I)) is called a minimax tra- 
jectory through x q . 

These are called the minimax elements of equations 

(1) and (3) with respect to x , y , and I. 

o o 

It may be desirable to define minlmax elements 
which do not depend on the time interval under consid- 
eration. Thus , supreme minimax elements are defined 
as follows. 

Consider the collection {1} of all closed intervals 

in I_, which contain t wherex(t ) =x andy(t ) =y . 
F o o o o o 

If there exist a control u* (t) e U and a point y* e R 

with the property that 

y* = max min max y ,i, u(t) ) , 

Ie{l} ueU tel 

then 

(1) y* is called a supreme minimax value of 
y(t,u(t}), 

(2) u* (t) is called a supreme minimax control, 
and 

(3) x (x , t; u* (t) ) is called a supreme minimax 
trajectory through x q . 

The two types of problems considered in this pa- 
per are as follows: 

Given equations (1) and (3) and sets U, X, and 

V 

(1) For t Q e I C Ij,, x q e X,y Q e R, andu'(t,I) 

e Ui show that u'(t,I) is a minimax control and find 
the other minimax elements. 

(2) For t e I„, x eX,y e R, and u* (t) e U, 

O f o o 

show that u* <t) is a supreme minimax control and find 
the other supreme minimax elements. 


hi view of these problems, several theorems were 
developed whichprovidesufficientconditions for mini- 
max control and supreme minimax control. The proofs 
of these theorems are presented in the next section. 
After that, the use of the theorems is illustrated by 
showing that certain elements are supreme minim x 
elements for a simplified booster control problem 
known as the constrained bending moment minimax 
drift problem. 

III. SOME GENERAL THEOREMS 

Some sufficient conditions for minimax control are 
presented in this section. Of the theorems presented 
here, the first two are probably the most important 
since the example in the next section depends heavily 
upon the use of these two theorems. Theorem 3. 1 
makes use of an integral cost function. To apply this 
theorem, it is important to express the cost as an in- 
tegral function, hence the expression (3) for the cost. 
An examination of (4) reveals that a control making 
y(t, u(t) ) attain its minimax value is equivalent to 
making 

t 

/ g(r, u(r) ) dr 

t 

o 

a minimax value since they differ by y , which is in- 
dependent of t and u(t). The four theorems are as 
follows: 

THEOREM 3. 1 : Given a real -valued piecewise 
continuous function g(t,u(t)) defined on a closed in- 
terval I and u(t) e U. If there exists u'(t) e U such 
that g(t, u'(t)) ^ g(t, u <t> ) for each t el ana each 
u(t) e U and if 

T 

|J g(t,u(t))dt| 
t 

o 

is bounded over I and U, then 

T 

min sup f g(t,u(t) )dt occurs for u(t) = u'(t) . 

u(t) eU Tel t 
' o 


Proof: Let u(t) be any control in U. From the 
typothesis there exists u'(t) el 1 such that g(t,u'(t) ) 
s g(t,u(t) ) for each t e [t ,r] CTI. Then it follows 
that 0 
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T T 

/ g(t,u'(t))dt s J g(t,u(t)) dt. 


T T 

Let / g(t,u(t)) dt-O(T), / g(t,u'(t)) dt-G'(r). 


Then G'(t) ^ G(r) for apy r s ij 

hence, sup G' (r) ^ sup G(t) for each r € I, 
Tel rel 

In terms of the original functions , 

T T 

sup f g(t,u'(t)) dtj£ sup f g(t,u(t) )dt 

rel t rel t 

o o 

for each u(t) e U, which is equivalent to 


max y(t,u(t)) = y(t ,u(t ))wheret e'C T (I ) • 

. . 'mm m Is 

t el 

because ylT.u^T)) is a minimax value for all ad- 
missible controls on I , It follows that 
s 

y(t ,u(t )) > max y(t,u(t)) because t does not 
m m t e i m 

s 

belong to I and max y(t,u(t)) s ’max y(t»ui(t)) 

8 tel tel 

s s 

because ui(t) is a minimax control on I . Further- 
more, by definition of u*(t) , 

max y(t,U!(t))= max y(t,u f (t)); 

tel tel 

s 

hence, max y(t,u(t)) > max y(t,u*(t)) which con- 
t el tel 

tradicts the assumption; hence, u f (t) is a minimax 
control on I. 

THEOREM 3,3: If l t C I 2 then y f (Ii) ^ y'(h)* 


T T 

sup f g(t,u f (t)) dt=min sup J g(t,u(t)) dt 

Tel t ueU t el t 

o o 

and the proof is complete. 

THEOREM 3.2 : Given a real -valued function 

y(t,u(t)) defined on a closed interval I and u(t) e U, 

If there exist a subinterval I C I and a control ui(t) 

defined on I such that there exists Tel where 
s s 

y(T,u t (T) ) = min max y(t, u(t)), 
ueU tel 

s 

and if there exists a control U 2 (t) defined on the com- 
plement of I in I (denoted C (I ) ) such that 
s is 

y(t, u 2 (t)) ^ y((T, ui (T) ) for eacht e Cjiy ,then 
there exists a minimax control on I given by 

fui (t) , t € I 
u'(t) =1 

[u 2 (t) , t e C^y . 

Proofs Assume u'(t) is not a minimax control on 
I. Then there exists a control u(t) c U onl such that 
max y(t,u(t) < max y(t,u*(t)). Furthermore, 
tel tel 


Proof: Let u'(t) be a minimax control on I 2 . As- 
sume that the maximum occurs at T£ I 2 . Then 

max y(t,u'(t)) = y(T,u’(T)). Therefore, 
tel 2 

y(T,u'(T)) s y(t,u'(t)) for eacht e I 2 . Sincel 2 I t , 
this expression holds for each t e I,, so, in particular, 

y(T,u'(T) ) > max y (t,u*(t) ) . It follows that 
telj 

max y(t,u'(t)) s max y (t,u f (t) ) or 
t e I 2 t e Ii 

y'(l2) — y'(Ii), and the proof is complete, 

THEOREM 3.4 : If u* (t) is a supreme minimax 
control with respect to x q and y Q whore xit^.u* (t > ) = 

x and y (t ,u*(t >) = y , then for any t * e I_, u* (t) 

O 0 O Oi? 

is also a supreme minimax control with respect to 
x(t Q ',u* (t ') ) and y(t o ',u* (t ’)). 

Proof: Let y'<I) be a minimax point with respect 
to I where I is any closed interval in I p containing t . 
Let y* be the supreme minimax value associated wifti 
u* (t). Since t Q and t ' are both in y, there exists 

I q € {i} which is sufficiently large to contain both; 

hence, y* => max y *(I) = max y'(I) = max y'(I') 
Ie{l} I5I o I*e{l»} 
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where the set of all closed intervals in 1^ con- 

taining t Hence, u* (t) is a supreme minimax con- 
trol with 0 respect to x(t ',u* (t ')) andy(t ',u* (t ')) 

for each t ' e I_. 0 0 

o F 


IV. APPLICATION TO A SIMPLIFIED BOOSTER 
CONTROL PROBLEM 


Application of the previous theorems to a booster 
control problem known as the constrained bending mo- 
ment minimax drift problem is presented in this sec - 
tion. This is the problem of minimizing the maximum 
lateral drift from a reference trajectory along flight 
time where the magnitude of the bending moment can- 
not exceed a critical value. 


be unrealistic since the object is to keep the maximum 
of |z(t)i as small as possible. Under these assump- 
tions , the problem under consideration is the following : 

Minimize |z(t)| subject to the conditions 

tei 

^ + C 1 <p + C 2 )3 = 0 

z - <K t *■ K 2 ) cp + K 3 /3 (7) 



Observe that the bending moment can, in this 
problem, be treated as a control variable in a very 
natural way. Define the new variable v in such a way 
that 


M v = M* <P+ M'B . (8) 

max ol R 

A. FORMULATION OF THE CONSTRAINED p 

BENDING MOMENT PROBLEM Using v as fhe control variable, the system (7) be- 

comes 


The rigid body equations of a typical vehicle 
are given in the following form: 

+ Cj ot + C 2 /3 = 0 

z = K t <p + K 2 a + K 3 /3 (5) 


where <p % a, /3, z represent angular deviation from 
reference , angle of attack, engine gimbal angle , and 
lateral positional deviation from reference, respec- 
tively. The wind contribution of the angle of attack is 

denoted by 01 and the vehicle speed by V. 
w 

The bending moment, given by a linear combina- 
tion of a and p 9 is constrained as follows: 

|M ! a + MJfil M fax 

I a fr\ max, (6) 

where M* and M* are bending moment coefficients and 

M is°ihe allowable bending moment. The problem 
max 

is to minimize the function 


<P = ix(p - Mv 

z = v<p + pMv (9) 

|v| * 1. 


where 


M = ( -Ci + C 2 


M' 


a 


M . ■>. ^ = (K 1 + K 2 

p 


K, 


M' 
o 

M* 

P 


P 


Ko 

— and M 


eg 

M^! M max* 

r 


The equations can now be puc in state variable 
form by making the substitutions x t = <p , x 2 = $ , and 

x« = ~ — . A constant K is defined by K = . The 

a v+p v+p 

constrained bending moment problem can thus be put 

in the following state variable form: 


Xj = x 2 

x 2 = pxi - Mv(t) ( 10) 

x 3 = Kx t + ( 1 — K) Mv(t) 


max 

tel 



v (t) | s i,osks1, 


subject to the subsidiary conditions (5) and (6). 

To obtain a model more amenable to analysis, two 
simplifying assumptions are made. It will be assumed 
that the wind angle of attack is zero and the quantify 

- will be neglected. The latter assumption may not 


where the object is to find a v(t) which minimizes 
max I x 3 (t) I . 

tel 1 1 

Since the purpose in this report is to illustrate 
the use of the theorems presented in the previous sec- 
tion on a minimax problem we will choose K=i. It is 
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felt that the same general techniques can be used on 
the more general problem (0 ^ K ^ 1) , but the analy- 
sis becomes more complicated. Choosing K = lyields 
the original version of the constrained bending moment 
minimax drift problem. In this problem the state var- 
iable equations are 


xi = x 2 


( 11 ) 


x 2 = fiXi - Mv(t) % n > 1, M>i (12) 

x 3 = xt, (13) 

where v(t) c V = {v(t) | | v (t) | ^ i} and each 

element in V is a piecewise continuous function of t. 
The quantities xj, x 2 , and x 3 represent attitude angle, 
its time rate of change, and lateral drift, respectively. 
The control variable v(t) represents a normalized 
bending moment. The problem is to find a control in 
V for which max | x 3 (t, v(t) ) j is minimized. The 
tel 

interval I is chosen sufficiently small such that pand 
M can be considered constant. 


Equations (11) , (12) , and (13) are in the form of 
(i) and (3) if we interpret y as | x 3 | and U as V. 

Let y = J x 3 1 = x 3 sgn x 3 for x 3 # 0, Then 
t 

y = x 3 sgn x 3 = x t sgn ( J xj dT + x 3 (t Q ) ) . 

t 

For this problem, equations (1) and (3) become 



* t 

y = x, sgn ( J Xidr + x 3 (t Q )) (15) 

*0 


where 

t 

f xi dT + x 3 (t ) * 0 
t 

o 

and (4) becomes 

t r 

y(t, v(t) ) = J x t sgn (f xjdp + x 3 (t ))dr +yU16) 

t t 

o o 


where xj is a function of v(t) through ( 14) . 

Consider in V the "extreme controls" v(t) = +1 
and v(t) = -1 on I. The trajectories x (t,x Q ; v(t)) of 
these two controls are seen in the phase plane as 


hyperbolas. The asymptotes of the hyperbolas inter- 
sect forming a closed region containing the origin 
( Fig. 1 ) . Let the region of admissible trajectories X 
be defined as the interior of this region. In other 
words , X is the region in which the control u = sgn x t 
provides stable trajectories. Formally, X is defined 
as follows: 

X = {x €R 2 | | Xj | + |^X 2 | < ■—} . 



THE REGION X 


FIGURE 1. THE TRAJECTORIES X(t,x ; ±1) SEEN 

in the f:iase PLANE ° 

B. SOME THEOREMS ON THE CONSTRAINED 
BENDING MOMENT MINIMAX DRIFT 
PROBLEM 

Some theorems dealing with the problem of 
findingminimaxvaluesof y(t, v (t) ) = |x 3 (t)| subject 
to equations (11), (12), and (13) are presented at 
this point. Let xj(t ), x 2 (t ) andx 3 (t 0 ) be denoted 
by x 10 , x 20 , and x 3 o, respectively. The first theorem 
is as follows: 

THEOREM 4.1: Let I* 10 ) e X and x 10 = 0. If 

1X 2 Q | 

x 20 * 0, then v(t) = sgn x< is a minimax control over 
any closed interval (t Q ,ti). If x 20 = 0, then v(t) = 0 
is a minimax control. 

Proof: 


(a) If x 20 = 0, then for the control v(t) = 0 
everywhere on Jt 0 , tj ] a solution to the differential 
equations ( 11) , ( 12) , and ( 13) is x t = 0, x 2 = 0, x 3 = 
x 30 on (t , t t ]. Hence, |x 3 (t)j =1x80118 certainly 
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a. minimal value of y(t, v(t)). Hence, v(t) = 0 is a 
minimax control under the condition x 10 = x 20 = 0, 

(b) Assume X 20 1 0, x 10 = 0 . 

Equations (11) and (12) can be inte- 
grated by variation of parameters to reveal the follow- 
ing solution for Xj: 

x 20 _ 

Xj (t, v) = sinh vjx (t -t ) 

•JJP 0 

M 4 

- f sinh«/J? (t-T) v(t) dr. (17) 

V7 t 

O 

i 

Letx 2 o > 0, For any t e [t ,t A ] , min Xj (t,v) oc- » 

0 veV 

curs for the control v(t) = +1, tc ft , t] . Thisisbe- 
cause sinh Vf?(t-T) 2 sinh«/^r(t-T)v(t) for v(t) =£ 1. 
For the control v(r) = +l, equation (17) becomes 

x 2 o _ 

x^tj+l) = — ~ sinh VjLi (t-t ) 

•J7 ° 

+ ~~ (l-cosh^/7i , (t-t o ) (18) 

Since! ^ I eX and x 23 > 0, if tj is sufficiently large 
\ x 20l 

there exists a unique T ^ t such' that x, ( T ,+i) = 0. 

P o P 

Furthermore, Xj(T ,+l) == xj(T ,v) for each v e V. 

P P 

For each t e [t ,T ) , the condition 0 ^ xj (t,+i) — 
Xj (t,v) holds. 0 It follows from equation (13) that 
x 3 (t,v) is monotonically increasing on ft Q , T^] for 
each v e V. It should be observed that for any 
T c f t Q , Tp ] themaximumof x 3 (t,+l) over the inter- 
val [t Q , T] occurs either at t=t n or at t=T. 



FIGURE 2. SKETCH OF x^t.+l) and x 3 (t,+l) 

Consider equation (16) over the interval ft 0 ,T] 
where y (t,v (t) ) = |x 3 (t,v(t))| . We will use The- 
orem 3. 1 to show that v(t) = +1 is a minimax control 
over [t^, T] where — T ^ T , We are given 


y =g(t,v(t)) where 

t 

g(t,v(t) ) = xj sgn ( J x t dT + x 3 (t Q ) ) where 
*0 

x t = xj (t, v) as g’ven by equation ( 17) . Equation ( 17) 
shows that xj (t,+l) ^ x 4 (t,v) where | v(t)| s i. 
t t 

For this reason, J x t (t,--1) dr == J Xj (t, /) dr 


from which it follows that 

t t 

sgn (/ x^r.+i) dr + x 3 (t ))<sgn( f x 1 (T,v)dr+x 3 (t )) . 

t t 0 

o o 

If the equality holds, then g(t',+i]t_ ^ g(t, r ) because 

x 4 (t,+i) ^ x t (t,v) . If the strict in qualify holds, 

then since (t,v) > 0 on ft , T] , it follows that 

g(t,--i) ^ 0, whereas g(t, /) > 0on(t ,T); hence, 

g(t,+i) ^ g(t,v). It has been established that, for 

t e [t Q ,T ] where t Q < T ^ T^, g(t, J-i) ^ g(t,v) for 

each v e V. It follows from Theorem 3. 1 that v(t) = 

+1 is a minimax control on ft ,Tl where Teft ,T 1. 

1 o’ 1 1 o p J 

The problem now is to find a minimax control 
over the entire interval ft ,tj] which reduces to 
v ( t) =+i on the subinterval ft ,T ]. Consider as a 
possibility the control v’(t) = sgn x£ For this control 
Xj given by equation (17) is a periodic function of t. 
In fact, x t (T + 1)= -x^t). Hence, y(t,/'(t)) is 
also periodic §nd never exceeds its maximum value 
on ft ,T ] ; that is, y (t,v ! (t) ) - max y(t,+l) foreach 
° P te[o,T ] 

P 

te(T , tj) . Hence, by Theorem 3.2, the control 
P 

v'(t) = sgn x 2 is a minimax control over the entire 
interval (t ,tj). A minimax control law has thus been 
found for x 10 = 0 , x 2 o > 0 which is independent of the 
final time tj. For the case x 10 = 0, x 20 < 0, an argu- 
ment similar to the preceding can be used. 

It can be shown further that , under the conditions 
of this theorem, the control v'(t) = sgn x 4 is, in fact, 
a supreme minimax control. This is stated as the 
following lemma. 

LEMMA: The control v = sgn xj is a supreme 

minimax control with respect to x e X and y e R 
, , „ . o o 

where x 20 * 0 and x 10 = 0. 

Proof: It has been established that v = sgn x 4 is 
a minimax control over rt 0 , t t ] for any real t,. Fur- 
thermore, the functions x(t, sgn xj) andy(t, sgn Xj) 
are both periodic. If the interval [t .tj] is chosen 
sufficiently large to include a complete period, then 
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'r ».**w,*wt* ,«* «•*»*%«* n 


max y(t, sgnxt) = max y(t, 3gnxj), 
te [t Q , ti] t € I 

where I is any interval in I T1 which contains [t ,tt] # 

if 0 

Hence 

max y(t, sgnxj) = max min max y(t,v) =y*. 
te[t o ,t!] Ie{l} veV tel 

The control v = sgn xjis thus a supreme minimax con- 
trol with respect to x q =( j 6 ^ x 30 e R where 

x 2 o *0. 

This lemma along with Theorem 3.4 establishes 
the key theorem for use on the constrained bending 
moment problem. 

THEOREM 4. 2 : The control v = sgn x* is a su- 
preme minimax control with respect to any point 
(x^ * 0) e X and y^ e R. 

Proof: Consider a point x e X such that x 4- 0. 
Consider the control v = sgn x^ and the corresponding 
trajectory x(t, sgn xj) through x . Let x(t^) = x , 

y(t ) =y . There exists t * el such that xWt *) = 0 
* o o of o 

and x 2 (t^ ! ) =£0, Hence, the previous lemma shows 

that v*= sgn Xj is a supreme minimax control with re- 
spect to x(t o *) andy(t^), It follows from Theorem 

3,4 that v = sgn x* is also a supreme minimax control 
with respect to and y^, 

C. SOME CONCLUSIONS ON THE CON- 
STRAINED BENDING MOMENT PROBLEM 

The final result is that v* = sgn xj is a su- 
preme minimax control for each nonzero x e X and 
x 3 (t ) e R. If the limit cycle trajectory is undesir- 
able, an approximate supreme minimax control can be 
obtained by setting v =sgn (x A + e x 2 ) where e is a suf- 
ficiently small positive number. 

The concept of supreme minimax control was in- 
troduced because of certain advantages it holds over 
an ordinary minimax control in the constrained bending 
moment problem. The principal advantages are its 
independence of the time interval I and of the initial 
conditions x q andy^. A supreme minimax control does 

not have the following shortcomings that some ordinary 
minimax controls may have: 

(1) A control which is minimax over I may, in 
general, lead to a poor control over some larger 


interval containing I, For this reason, piecing to- 
gether minimax controls over various intervals may 
not lead to a minimax control over the union of these, 

(2) Minimax control over a fixed I can lead to 
different control laws for different initial conditions x 
andx 3 (t o ), 

(3) For one set of initial conditions, a minimax 
value can be less than for another admissible set of 
initial conditions. The most practical interest, how- 
ever, is in the "worst" minimax value. 


V. A SUMMARY OF SEARCH PROCEDURES FOR 
MINIMAX CONTROL 


An outline of methods for seeking a minimax con- 
trol or supreme minimax control which are based 
primarily on Theorems 3. 1 and 3. 2 is presented for 
a class of minimax problems. 

No claim is made that these techniques can be ap- 
plied in each case. They are presented merely as 
suggestions as to |how Theorems 3. i -3.4 might be 
used in constructing a minimax control law. 

Consider the class of problems in which equation 
( 1 ) is linear in x and u. 

x = Ax+bu(t) +Cw(t), x(y =x q (19) 

y = h(t,x,u(t) ) , y(t Q ) =y Q (20) 

In this problem, w(t) is a known function of time, A 
is an nxn matrix, b and c are n-vectors. Given a re- 
gion X and x for which there exists at least one con- 
trol in U such that the solution trajectory remains in 
X for all tel and given y Q e R, the problem is to find 
a minimax control over 1. A possible approach is as 
follows: 

(1) Solve equation (19) by the method of varia- 
tion of parameters and obtain x = <p(t,u(t) ) . For this 
solution equation (20) becomes 

= h(t,<p(t,u(t)), u(t) ) =g(t,u(t)), y(t Q ) =y Q (21) 

( 2) Pick a subinterval I C I in which the peak 
value of y is expected to occur. For example, if •« (t) 
has a peak, for many physically meaningful cost func- 
tions , the peak of y might be expected to occur near 
the peak w(t). In these cases, I should contain the 
peak time of w(t). 
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VI. CONCLUSION 


(3) Find, if possible, a control uj(t) € U such 

that g(t,ui(t)) ^ g(t,u(t) ) for any t e I andu(t)eU. 

Then by Theorem 3.1, ui(t) is a minimax control on I . 

Lety’(I ) denote the minimax value on I . 

^ s s 

(4) Find any control U 2 (t) e U defined on Cj(I g ) 

for which the peak value of y ( t, u 2 (t) ) is less than or 

equal y'(I ) . Then by Theorem 3.2, the control 
s 

{ ui (t) , t e I g 
u 2 (t) , t e Cjiy 

is a minimax control over I. 

(5) If possible, pick ug(t) such that the function 
u*(t) is independent of the interval I. Inthiscase, u’(t) 
is a supreme minimax control. 


Control problems in which the object is to mini- 
mize the peak value of a functional constitute a very 
practical class of problems having many applications 
in aerospace vehicle guidance and control. In this re- 
port some sufficient conditions for solution of problems 
in this class were stated and applied to a simplified 
version of a booster control problem. Some conjec- 
tures for application to other problems were also pre- 
sented, It is hoped that further study in this area will 
lead to more systematic procedures for determining 
minimax control elements. 
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V. GUIDANCE STUDIES 



CONCEPTS OF THE ITERATIVE GUIDANCE LAW FOR SATURN LAUNCH VEHICLES 

By 

Clyde D. Baker 



SUMMARY 

- / ; ' 


| Beginning with the V2, rockets in the mid- 1940*5 
and continuing through the Redstone, Jupiter and 
Pershing missiles, the delta minimum concept was 
employed successfully. The delta minimum conc?pt 
requires that both nominal and perturbed trajectories 
fly essentially the same geometrical trajectory regard- 
less of other considerations. Analog computation was 
employed for the execution of the delta minimum con- 
cept. 


In 1960, research work was begun to develop new 
guidance concepts for the Saturn space vehicles. This 
work was motivated by the development of new math- 
ematical techniques for maximization of payload 
through optimization methods and by the development 
of digital computers to replace analog computers as 
on-board hardware. It was also obvious that space 
trajectories would require greater flexibility to cope 
with sudden changes, such as engine out, and that more 
flexibility must be permitted in the selection of flight 
profiles. 


The accomplishment of any of these missions 
requires proper functioning of the control, navigation, 
and guidance systems. Since these functions are inter- 
coupled, in the final analysis they must be considered 
simultaneously. However, for the purpose of eApla- 
nation we will discuss each as a separate function. 

The control system then has as its primary func- 
tion the supervision of the rotational motion of the 
vehicle, subject to the constraint of the prescribed 
reference attitude while maintaining structural loads 
within acceptable bouids in the presence of disturb- 
ances such as wind. Pilch and yaw attitude angles are 
prescribedby the guidance system and are then intro- 
duced as reference angles into the control system. 
The control system then rotates the vehicle into the 
prescribed reference direction hile simultaneously 
insuring that maximum allowable bending moments are 
not exceeded and that propellant sloshing does not 
reach disastrous amplitudes. Simply stated, the 
control system maintains the desired vehicle attitude 
angles and angular rates while preventing the struc- 
tural destruction of the vehicle. 


This paper presents the basic concepts of the 
iterative guidance law developed for Saturn launch 
vehicles to meet these new requirements of space age 
guidance. The development of the iterative guidance 
law and the results and ideas presented in this paper 
are due primarily to Mr. Helmut J. Horn and his 
associates in the Dynamics Analysis and Flight 
Mechanics Division of the Aero-Astrodynamics 
Laboratory, 


The purpose of the navigation system is to provide 
accurate information describing the current position, 
velocity and acceleration of the vehicle as a function 
of flight time, or some other convenient variable. 
This information is provided with respect to some 
reference coordinate item. Thus it may be position, 
velocity, and acceleration of the vehicle with respect 
to the launch site or it may be given with respect to 
the target point. Such information may be provided 
by an inertial platform together with integrating accel- 
erometer as in the Saturn vehicles or it may be provided 
by radar measurements, or visual tracking, etc. 


I. INTRODUCTION 


A typical space guidance problem is that of placing 
a space vehicle into a specified circular orbit around 
the earth. The particular iterative guidance solution 
of this problem is presented in this paper to illustrate 
the basic principles and concepts involved. Among 
many other applications, iterative guidance will inject 
a spacecraft into a specified lunar orbit or will soft 
land a vehicle at a specified point on the surface of the 
moon. The discussion here will be limited primarily 
to the problem of injection into a circular orbit about 
the earth since the principles of the guidance scheme 
are essentially the same for all specific applications. 


Finally, it is the purpose of the guidance system 
to use the position, velocity, and acceleration infor- 
mation provided by the navigation system to determine 
the reference attitude for the control system. This 
means essentially that, on the basis of known position, 
velocity and acceleration of the vehicle, the guidance 
system prescribes the required direction of the thrust 
vector necessary to achieve the mission of the flight. 

In the beginning of the Saturn guidance studies, it 
was assumed that guidance would become active only 
after the burnout of the first stage. Trajectory shaping 
of the first stage was constrained primarily to prevent 
excessive aerodynamic structural loads and aero- 
dynamic heating. Within these constraints , a zero-lift 
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first stage, trajectory was then computed to maximize 
tiie payload when connected to an optimized second 
stage toorbit. More recently, angles of attack during 
the first stage have been permitted in order to further 
increase the payload delivered to orbit. 

For the specific guidance problem of entering a 
circular earth orbit discussed here, it is necessary to 
cut off the engines when the specified altitude is reached 
with zero degrees path angle and circular velocity. All 
three of these conditions must be satisfied simultane- 
ously. 

To illustrate the solution of the iterative guidance 
problem, the complete equations of motion of the 
vehicle are presented first. These equations are then 
simplified step by step to illustrate the essential ele- 
ments of the problem. 


Since guidance is initiated after the vehicle is out 
N A 

of the atmosphere, — and — are both zero and the 
m m 

equations may be simplified to 


x = f/ 1, ~ + gj 

gl = Gj^g 

(4) 

y = F 2 (1, m + g2 

g2 = Ga. (1) g 

(5) 

*"F 3 (l, m + 6 

gs = G 3 (l) g . 

(6) 


The values of the thrust force unit vector compo- 
nents are given by 


II. EQUATIONS OF MOTION 


The equations of motion of the Saturn vehicle are 
usually written in the following form: 


F^ = cosx cost cos/3^ cos/3 ^ - sinx cos /3^ sin/3^ 


- cosx sin/3^ sinr 


(7) 


x = F^- + n} 1 — - C 1 < 1) - + GjW g 
1 m m 1 m 

y = F 2 (l) — + N 2 ( ^ — - C 2 (l) — + G 2 (l) g 


z = F 3 (i) ^ + n 3 (1) ^ - c 3 (1) ^ + g 3 (1) g 


( 2 ) 


(3) 


^ ^ F 2 ^ = sinx cost cos/3^ cos/3^ - sinx cos/3 ^ sin/3^ 


+-cosx cos/3^ sinT 


( 8 ) 


F 3 ^ = sinT cos/3 cos/3 + cost sin/3 (9) 

| g b 


where 


where 

F^ is aunitvectorinthe direction of thrust force, F, 
N^ is a unit vector in the direction of normal force, N, 
d 1 * is a unit vector in the direction of axial force, A, 


X is the vehicle attitude angle in the pitch plane, . 
t is the vehicle attitude angle in the yaw plane, 

/3^ is the engine gimbal angle in the pitch plane, and 
/3^ is the engine gimbal angle in the yaw plane. 


is a unit vector in the direction of gravitational 
force, g, and 

rn is the mass of the vehicle. 

Subscripts l, vs, and 3 indicate components of the 
unit vectors along the x, y, and z coordinates, respec- 
tively. These coordinates have their origin at the 
center of the earth with the y-axis parallel to the 
gravity vector at launch. The x-axis is parallel to the 
downrange direction and the z-axis completes a right- 
handed system. 


Equations (7), (8), and (9) are exact represen- 
tations of the r espec tive dir ec tions of thrust. However , 
since these equations are coupled in the pitch and yaw 
axes, the above definitions of x. t, ft , and /3„ are not 

strictly correct. Because of the coupling, they may 
not correctly be referred to as being in the pitch and 
yaw planes. Using the fact that the angles t, /3. and 
/3^ are small, the definitions given above are approx- 
imately true and allow a visualization of these equa- 
tions. 


t 


57 





PAGE 


* 


The values of gj, g 2 , and g 3 are somewhat compli- 
cated functions of the vehicle altitude and latitude and 
will not be given explicitly here. 

Equations ( 4) , ( 5) , and (6) are the precise equa- 
tions of motion of a vehicle in vacuum on a rotating 
spheroidal earth and may include a full six degrees of 
freedom with the inclusion of the dynamics equations 
relating gimbal angles and attitude angles. These 
equations may also be simplified at this point to three 
degrees of freedom by making certain simplifications 
concerning the dynamics of rotational motion. For 
guidance studies, the latter choice is adequate. 

Accordingly, the dynamics of angular motion will 
be neglected, and the thrust direction will be assumed 
along the longitudinal axis of the vehicle. Further- 
more, guidance will be assumed to restrict the 
trajectory to a plane so that the equations of motion 


may be reduced to 


.. F 

x = — cosx(t) + g! 

(10) 

y = m sinx(t) + g2 • 

(11) 

For injection into a circular orbit about the earth, 
the following conditions must be satisfied. 

x 2 +y 2 = r 2 
c c c 

(12) 

x 2 +y 2 = v 2 
c c c 

(13) 

. • 

xx + y y =0 
c c J c J c 

(14) 


where the subscript c refers to engine cutoff time at 
injection into orbit. 

The guidance problem then is to prescribe the 
angle x(t) so that equations (12), (13), and (14) are 
satisfied simultaneously at injection. This problem 
as stated does not have a unique solution. There are, 
in general, many trajectories and corresponding 
steering angles x(t) which will fulfill these end condi- 
tions given by equations ( 12) , ( 13) , and ( 14) . 

By imposing the requirement that these terminal 
coj 'itions be fulfilled in minimum time, or equiva- 
lently, that they be fulfilled by minimum fuel con- 
sumption, the steering law becom es unique, and may 


be solved by the use of calculus of variations. How- 
ever , the calculus of variations solution to this proble'm 
is too complicated to be carried out in flight by an on- 
board computer. Consequently, further simplifica- 
tions are made to reduce the complexity of the problem 
while maintaining the desired accuracy and minimum 
fuel consumption. 


ID. THE NATURE OF THE MODIFICATION 
OR SIMPLIFICATION 


As pointed out earlier, the calculus of variation 
solution of the optimum thrust direction to maximize 
payload on a rotating, spheroidal earth requires com- 
putation facilities far in excess of those available with 
an on-board computer. 

The iterative guidance concept derives the optimal 
thrust direction law for a flat earth model with constant 
gravitational field and then applies this law in the 
realistic earth model. It is important to note that 
only the optimal control law from the simplified model 
is used. The equations of motion for Saturn are 
calculated on a rotating spheroidal earth. 

Obviously, a careful check must be made to deter- 
mine if the optimal guidance law derived for the 
simplified model provides optimal performan . and 
accuracy for the actual Saturn trajectories. Many 
numerical calculations demonstrate that the optimality 
can be preserved. 


IV. THE ITERATIVE GUIDANCE MODE CONCEPT 


Since equations ( 10) and (11) are too complicated 
for on-board computation, further simplifyingassump- 
tions are made that 


gi = 0 (15) 


g = constant ( 16) 

and equations (10) and (11) then become 

x = cosx ( 17) 

m 

•• F 

y =~ sinx - g. ( 18) 
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Equations (17) and (18) may be written as a 
system of first order equations in the form 


x A = ^cosx = x 

x 2 = Xj = x 
• F . 

x 3 = “ sinx - g = y 


x 4 -x 3 


= y 


(19) 

( 20 ) 

( 21 ) 
( 22 ) 


To apply the maximal principle to this system of 
equations, define 


x 5 = 1 


where the quantity to be minimized now is 


(23) 


J - J x 5 dt , 


subject to the end conditions 


x(t ) unspecified 

x 2 (t c ) unspecified 

(24) 

y(t c> = y te 

X * ( V = X *tc 

(25) 

i(t c )=x tc 

x i ( V =x ‘tc 

(26) 

y(t c )=y tc 

x 3 ( t c )=x 3tc t 

(27) 


These equations are the flat earth equivalents of 
equations ( 12) , ( 13) , and ( 14) . 


The Hamiltonian of this system of equations then 


is 


F F 

H = p i “ cosx + p 2 x, + P 3 ( — sinx - g) 


+ F 4X3 - 1 


(28) 


To maximize the Hamiltonian wUh respect to the 
control variable x» calculate 


dH F F 

— = - P x — sinx + P 3 cosx = 0 
9X 4 in 0 m 

which implies 


tanx = 


From the canonical equations, 

• -0H 

Pl " 3x t “ ~ P * 


p ’=f=° 


Pl 3x,- p < 


• -0H 

P -^=° * 


Hence, 

P 1 = _P 20t + P 10 


P 2 ~ P 20 


P 3 _ -P 40t + P 30 

P4 = P40 • 

Thus , equation ( 30) may be written 

-Pii't + Pan 


tan* = 


-P 2 ot + p io 
From the transversality conditions 


( 2 Pi ) dx =0 

i = i _ i 


x i =x i° 


Since x(t c ) is unspecified, P 2 o = 0 
( 39) may be written 


(29) 

(30) 

(31) 

(32) 

( 33) 

(34) 

(35) 

(36) 

(37) 

(38) 

(39) 

(40) 

and equation 
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tanx = ~^ 4l>t p P|ift = A + Bt 
p lo 


(41) 


where 

A=| 311 (42) 

P 10 


B=^ 4(l (43) 

P 10 

Equation (41) then is the optimum steering law to 
maximize payload for a trajectory on a flat nonrotating 
earth with constant gravitational field. This expres- 
sion is further simplified by assuming that 


X = a + bt . 


(44) 


This expression may be substituted into equations 
(17) and (18) to obtain 


• • 

x = 


F 

— cos (a + bt) 
m 


(45) 


Equations ( 47) , (48) , and ( 50) may be substituted 
into equations (26), (27), and (25), respectively, to 
obtain three simultaneous equations in the three 
unknowns a, b, and t c . 

The simultaneous solution of these three equations 
provides the values of a, b, and t c required for cal- 
culating the optimal payload for the flat earth case. 
This is an open loop solution to the problem. This 
process may be made a closed loop process by re- 
peating it and recomputing new values of a, b, and t Qi 
at each second or so of flight. This is the desired 
solution sought from the beginning of the problem. 
Using the values of a, b, and t c , the resulting value 
of x is taken as an approximate solution to equations 
(4) and (5). This is the basic concept of the iterative 
guidance mode. The remainder of this paper is con- 
cerned with three topics: 

(a) Methods for solving the three simultaneous 
transcendental equations in three unknowns. 

(b) Methods of applying the flat earth optimal 
solution to the realistic case and improving the opti- 
mality in this application. It will later be shown that 
accuracy is not a problem because of the sequence of 
computation. 


y = 


F 

— sin ( a + bt) 
m 


(46) 


(c) Stability and accuracy of the iterative guid- 
ance mode. 


Integrating equations ( 45) and ( 46) yields 


. r c F 
x = I — cos ( a + bt) dt + x 
c J m o 


and 



F 

— sin (a + bt) dt - gt + 
m 



(47) 

(48) 


V. SOLVING THE SIMULTANEOUS EQUATIONS 


Returning to the transversality condition of equa- 
tion (40), assume for the moment that y(t c ) is 
unspecified as well as x(t c ). Then both P 20 and P 40 
are equal to zero, and the optimum steering law is 
given by 


and integrating equations (47) and (48) yields 


tanx = 


Pan 

PiO 


= X a constant. 


(51) 



( 49) Equation ( 51) states that, for aflat earth, the optimum 
steering program to achieve a given velocity, without 
position constraints, is a constant. 


and 


For this case, equations (45) and (46) can be 
written 



sin (a + 


bt) dt 2 - 

La 


+ y 0 t+y 


0 


(50) 


• * F rv 
x = — cosx 
m 


(52) 
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y = — sinx - g • 


Integrating equations ( 52) and ( 53) 


x = v in cosx + x 

c ex T - T o 

o cl 


y c*v'“ sln * +y cr g ‘c 

\ o cl 


where 


(53) 


(54) 


(55) 


ex 

(56) 

T - t 

mill 

~ ft 

(57) 

m 


= I g = characteristic velocity. 

(58) 


Equations ( 54) and ( 55) may be solved for x 


X = arc tan 


• ft 


y c " y o +gt c 


x - x 
c o 


(59) 


It may be observed that this same argument holds 
without the introduction of % in the first place. The 
introduction of this angle in the manner described 
makes possible the further simplification described 
below in equations ( 60) , (6i) , (62) , and ( 63) . 

Equations (45) and (46) can now be written as 


• • F ~ 
x = — cos ( X - kj + k 2 t) 


y = ~ sin ( x - kj + k 2 t) 


('30) 


( 61 ) 


Simplifying, assuming that ( -kj + k 2 t) is a small 
angle , 


cosx ** cosx + kj sinx - k 2 t sinx 
sinx ** sinx - k t cosx + k 2 t cosx 


(62) 

(63) 


Substitute these values into equation (61) and 
integrate to obtain 


y c - ^i(ki» k 2 , t^) , 


Integrate once again to obtain 


(64) 


Equation (59) gives the value of the constant 
steering angle x which is required to guide a vehicle 
from its initial point to specified values of velocity 
x(t c ) andy(t c ) in minimum time. These are the same 
values of velocity components which are required when 
the position constraint is imposed where, of course, 
the time will be increased. This angle x represents 
the major part of the angle x so that x may be divided 
into its principal part x Plus a second smaller part 
(-k 4 + k 2 t). 

By adding this angle (-k t + k 2 t) , x is changed so 
that, in principle, it might not be expected that the 
same terminal velocity components are obtained after 
the addition. However, since two parameters are 
added, the two may be chosen simultaneously so that 
the velocity requirements may be preserved while 
meeting the position constraint. 


y = f 2(ki» ka, t ) . (65) 

c o 

Equate (64) with (55) and solve the resulting 
equation simultaneously with equation ( 65) for kj and 
k 2 . 

Th’is, the equation x = 5? “ ^ + k a t, where x » kj, 
and k 2 are derived in this manner, fulfills both alti- 
tude and velocity conditions. 

To solve for t c , consider the velocity to be gained 
at a given instant of time. 

Av = [(x o - x) 2 + (y c - y +gt c ) 2 J (66) 

together with the characteristic velocity equation 
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t =T[l- e - av/v ex] . 
c 


(67) 


Equations ( 66 ) and (67) may be solved simulta- 
neously for t c . Thus, kj, k 2 , and t c have been obtained 
for the flat earth case. 


VI. GUIDANCE OVER A SPHERICAL EARTH 


The methods employed in the derivation of the 
guidance equations for a spherical earth are essentially 
the same as those used with a flat earth model. How- 
ever, the equations of motion must be modified to 
account for the constantly changing magnitude and 
direction of the gravity force. An average gravity 
magnitude g* , and an average gravity direction 0 * , 
between the current point on the trajectory and the 
final injection point are updated at the beginning of 
each computation cycle. 


For deriving the guidance equations on a spherical 
earth, it is convenient to change coordinate systems 
so that the terminal conditions expressed in equations 
(12), (13), and (14) may be expressed in a more 
convenientform. The new coordinate system is chosen 
so that the r\ axis coincides with the local vertical at 
the injection point. The £ axis is chosen parallel to 
the local horizontal at injection. In this coordinate 
system, the equations of motion are 


.. F 

t = — cosx „ - g* sin 0 * 
s m £ 

( 68 ) 

.. f 

ti = — sinx + g* cos 0 * . 

' m £ 

(69) 


The terminal conditions in this coordinate system 


are 


£ c unspecified 

(70) 

| = v 
s c c 

(71) 


(72) 

d 3 

11 

a 

(73) 


These equations may now be solved exactly as in 
the x, y coordinate system. 


VII. SUMMARY OF GUIDANCE EQUATIONS AND 
REQUIRED INPUT CONSTANTS 


For the calculation of the guidance of a single 
stage to earth orbit, the input constants are 


• • 


V ex’ V V k andt c ' 
c c c 


The initial measured values are 


— , x , y , x , y , v , 0 
m’ o J o o J o o Y o 


Calculate 


T = V /— 
ex m 


A < = -v in (1 - — ) 


<h 


ex 


Ao = A 1 t - v t 


ex c 


= -Ao + t A 4 
6 c 1 


V t 


= Ao T 


ex c 


■ 2 
A 3 + v_t 


o c 




f t 


= x cos (0 + 0 , ) - y sin (0 + d> ) 
o o t o o t 

c c 


= x sin ( 0 + 0 ^ ) + y cos (0 + 0 X ) 
o o t o o t 

c c 


L = K cos + ) - y„ sin <4* + $ 4 . ) 

^o o o t o o t 


\ = 5 i o Sin(0 o + <#, t > + 3 ? o COS ( 0 o +<?> t > 


g* = i ( s 0 + s t ) 


4>* - h 
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Av* = (i t - l 0 - g* t c sia ^ ) 2 

V ° 

+ (m -« + g*t cos**) 2 ^ 
't 'o o 
c 

Av a = -v n ( 1 - t /t) 

0 ex o J 


Solve simul- 
taneously for t 


for stability and accuracy. The F/m derivative is 
small during the entire flight, eliminating this usually 
rather noisy measurement as a trouble source. How- 
ever, because the trajectory optimization is based on 
a predicted relation of the future thrust profile for a 
stage to the instantaneously measured value, any 
major thrust change will cause a performance loss. 



= tan 1 


r? t - T) 0 + g* t c cos <j>* 

T t - l 0 -s* ‘o 8 ^* 

c 


P = A 3 coax 

Q = A 4 cosx 

R = \ * \ ~ Vc + ^ l c 2 cos< * > * " A s sin X 

c 


ki 


— AaR. ... 
AjQ - A 2 P 


k , = M_ 

2 AjQ - A 2 P 


The other derivatives start at low values and 
increase approximately inversely proportional to the 
time-to-go (for velocity errors) or its square (for 
displacement). The tightening of the guidance loop 
toward the end of flight is very desirable because it 
keeps residual errors small. However, it creates a 
potential stability problem. This problem was elim- 
inated without causing a significant error by stopping 
computation of the steering equations at a given time- 
to-go (e. g. , t c = 20 seconds) and flying open luop. A 
better method is to freeze the time-to-go at a minimum 
value and continue guidance. 

i 

The low guidance gains at early flight make the 
system very tolerant to major disturbances, noise and 
time lags during this phase. 

Guidance scheme errors for realistic variations 
of initial conditions (Table I) are very small. The 
effects of performance variations, changes in air 
density, and winds are equally insignificant. 


X £ = * I “ (k l " k 2 4) 



This angle x is introduced into the control system, 
and the computation cycle is started over again. As 
the computation progresses through the flight, com- 
putation of kj and k 2 is stopped when t c < 20 seconds. 
Cutoff is given when t<. = 0. The system of equations 
listed above is calculated repeatedly throughout the 
flight, treating each time point as new initial condi- 
tions. Errors of the approximation introduced tend to 
zero at terminal time, and terminal conditions are 
fulfilled precisely. 


VIII. STABILITY, ACCURACY, AND OPTIMALITY 
OF THE I. G. M. 


TABLE I. GUIDANCE SCHEME ERRORS FOR 
REALISTIC VARIATIONS OF 
INITIAL CONDITIONS 


[initial Stage Variable 

Payload 

Loss 

percent 

Injection Errors 

AXj 

km 

Axt 

m/s 

Ayi 

km 

Ajh 

m/s 

Altitude 

m 

Velocity 

m/s 

Path Angle 
degrees 

+2.7 

0 

0 

0 

.11 

,1 

0 

+.001 

0 

+143 

0 

0 

.32 

.13 

04 

0 

0 

- 57 

0 

(► 

.07 

.1 

+.01 

+.001 

0 

0 

1.0 

0 

.11 

.1 

0 

+.001 

0 

0 

0 

+78 

.14 

.1 

-. 01 

0 

0 

0 


-80 

.11 

.1 

0 

+.001 


A time lag of 5 seconds from measurement to 
steering command causes no error and no loss of 
weight in orbit. A 40 second lag caused 3 km altitude 
error and li percent payload load loss. 


The partial derivatives of attitude with respect to 
the state variables are the most significant criteria 


Periodic thrust fluctuations with a maximum 
amplitude of 65 percent of nominal and periods of 5 - 
100 seconds create no serious stability problem. 
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DEVELOPMENT QF A LASER DOPPLER FLOWMETER FOR GAS VELOCITY MEASUREMENT 

By 

W. Foreman- and Robert M. Huffaker 


that application of ordinary analytical techniques to the 
reverse flow phenomenoncannot be made with complete 
confidence. Experimental measurements of the local 
temperature, density, and flow velocity of the exhaust 
gases in the vicinity of the base region are required 
to supplement the analytical treatment. 

Reverse flow phenomena also occur in cold flow 
models of rocket engines — that is, scale models of 
rocket nozzle clusters in which ordinary air at room 
temperature is forced through the nozzles at high 
speed. Both the analytical and experimental ap- 
proaches to the reverse flow problem are simpler in 
the cold flow case. For this reason, it is logical to 
begin a study of reverse flow phenomena with a cold 
flow model. 

This paper is concerned with the development of 
anew experimental technique for measuring one of the 
variables associated with reverse flow fields in cold 
flow models, namely, local flow velocities in the 
neighborhood of the base region. 


II. THE LASER DOPPLER FLOWMETER 


The following are the two principal conditions 
which a cold flow velocity measurement must satisfy: 

( i) It is necessary to measure local flow veloci- 
ties with treasonable accuracy. 


SUMMARY 


The measurement of gas velocity for a gas that 
has variations in pressure and density has proven to 
be very difficult using the more conventional methods 
of probes such as electron beams, gamma rays or 
x-rays. These methods require the detecting appa- 
ratus to be exposed to the radiation for long periods 
of time, or they have pressure limitations. The large 
magnitudes of the Doppler shifts obtainable, using a 
CW gas laser, suggest that local velocity vectors can 
be measured by using the visible light from the laser. 
Optical heterodyning of the laser Mght scattered from 
the flowing gas with a portion of the incident laser 
beam split off from a beam splitter will produce a beat 
signal at the frequency of the Doppler shift due to the 
motion of the gas. Measurement of the beat signal 
frequency , together with the geometry of the optical 
system, determines the flow velocity of the gas. A 
probeless velocity measuring technique has been de- 
veloped based on this principle. This technique has 
been checked in a wind tunnel at speeds up to 500 
ft/sec (150 m/sec). Work is continuing to determine 
velocity limitations. This technique also gives a dis- 
tribution of velocities due to the turbulence level of the 
gas flow. This technique is also applicable for the 
measurement of velocity of any moving surface. Some 
applications will be the measurement of panel flutter 
and for a calibration of dynamic pressure transducers 
or other applications where knowledge of particle flow 
is desired, f 
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I. INTRODUCTION 


(2) The pei'turbation induced in the flow pattern 
by the measurement technique must be negligible. 


Overlapping of the exhaust plumes from the vari- 
ous exhaust nozzles of high-thrust rocket engines gives 
rise to an undesirable flow characteristic called "re- 
verse flow" in which some of the hot exhaust gases 
flow back toward the base of the rocket. Impingement 
of the hot reverse flow on the base region causes local 
heating. The effect becomes more pronounced at high 
altitude, where the exhaust plumes expand due to the 
low ambient pressure. To achieve a satisfactory 
thermal balance in the base region of the rocket, the 
heating effects due to reverse flow must be taken into 
account. Unfortunately, the problem is so complex 


In addition to these requirements, it is highly desirable 
to measure the complete velocity vector v* at a given 
point, rather than a single component of v* along a 
given direction. This could be done, for example, by 
measuring simultaneously the components of v* along 
three noncoplanar axes which intersect at the given 
point. 

None of the gas flow velocity measurement tech- 
niques currently in use fully satisfy the requirements 
listed above. One promising approach to the problem 
which has rotyet been seriously exploited is the optical 


* Employed by Brown Engineering Co. , Huntsville, Alabama 


60 


Doppler technique. In this approach, use is made of 
the fact that when monochromatic light is scattered 
from moving particles, the frequency of the scattered 
light differs from that of the incident light (Doppler 
effect) . The Doppler frequency shift depends both on 
the velocity of the moving particles and on the geom- 
etry of the scattering. If the scattering geometry is 
fixed, measurement of the Doppler shift gives suffi- 
cient information to determine the velocity of the mov- 
ing particles. 

Before the development of the continuous -wave gas 
laser, there was no monochromatic light source with 
sufficient intensity to serve as the source in an optical 
Doppler flowmeter. Furthermore, even if a suitable 
source had been available, the best spectrographic 
techniques could not have furnished the resolution 
necessary to measure the small Doppler shifts asso- 
ciated with ordinary flow velocities. For example, if 
monochromatic light of wavelength \ =6000 A is re- 
flected at normal incidence from a body moving at a 
speed of 10 ft/sec (3 m/sec) the resulting Doppler 
shift is 10 MHz, which corresponds to a wavelength 
shift of approximately 0. 0001 A. A wavelength shift 
this small is well beyond the range of ordinary spec- 
tographic resolution. 

The development of the continuous-wave gas laser 
furnished a light source which meets the stringent re- 
quirements of the optical Doppler flowmeter. The 
concomitant development of the technique of optical 
heterodyning furnisheda Doppler shift detection mech- 
anism which possesses the necessary resolution. In 
optical heterodyning, two coherent, monochromatic 
light beams with different frequencies are directed 
simultaneously onto the photosensitive surface of a 
photodetector (for example, the photocathode of a 
photomultiplier tube) . The intrinsic square -law re- 
sponse of the photodetector leads to an electrical out- 
put signal the frequency of which is equal to the fre- 
quency difference between the two incident light beams . 
Optical heterodyning thus allows a direct measurement 
of the frequency difference between two coherent, 
monochromatic lightbeams from a few Hertz up to the 
maximum frequency response of the photodetector 
(typically several hundred mega -Hertz for fast photo- 
multiplier tubes) . The development of these two basic 
tools brought the optical Doppler flowmeter into the 
realm of feasibility. 

A laser Doppler flowmeter designed to measure 
liquid flow rates has been operated by Yeh and 
Cummins [1] at the Radiation Laboratory of Columbia 
University, A simplified diagram of the Yeh-Cummins 
laser flowmeter is shown in Figure i. Polystyrene 
spheres of diameter 0. 557 n were placed in water to 
form a dilute colloidal suspension (1 part solid to 
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FIGURE 1. A SIMPLIFIED DIAGRAM OF THE 
LASER FLOWMETER OF YEH AND CUMMINS 


30,000 parts water by volume) • These spheres served 
as optical scattering centers for the light beam gen- 
erated by a He-Negajs laser operating at a wavelength 
of 6328 A. The continuous power output of the laser 
was approximately 5 mW in the TEM 0 o mode. The 
contaminated water flowed from a large, constant-head 
reservoir through a flow tube of circular cross section. 
The flow rate was controlled by an outlet needle valve. 
The laser beam passed axially down the flow tube, and 
the light scattered at an angle 0 by the scattering cen- 
ters lying within a small linear segment of the laser 
beam was collected by an optical system, combined 
with a fraction of the incident laser beam by a system 
of beam splitters, and focused on the photocathode of 
a photomultiplier tube. The scattered light suffered 
a Doppler shift of frequency given by 



2 v n 
A 

o 



(i) 


where n is the index of refraction of the liquid, v is 
the flow velocity, and \ = 6328 Ais the vacuum wave- 
length of the incident laser beam. Optical heterodyn- 
ing of the scattered and unscattered beams on the pho- 
tocathode gave rise to a sinusoidal signal in the 
photomultiplier output at the difference frequency f^. 
Measurement of this frequency allowed the flow veloc- 
ity v to be calculated from equation (1). By translat- 
ing the beam back and forth across the diameter of the 
flow tube, Yeh and Cummins were able to measure the 
velocity profile across the tube. They obtained a 
parabolic "elocity distribution in excellent agreement 
with the theory of laminar flow. The minimum detect- 
able velocity in this experiment was estimated to be 
about 4 x 10“ 3 cm/sec. 


An experimental feasibility study was initiated 
within the Brown Research Laboratories to verify the 
feasibility of applying the laser Doppler flowmeter to 
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gas flow fields. In May 1965, the Doppler' shift due to 
scattering of visible gas laser radiation from a con- 
taminated flowing gas was detected for the first time 
by Foreman, George, and Lewis [2], thereby verifying 
the basic feasibility of the laser flowmeter. Further 
experiments have since been conducted to verify the 
feasibility of applying the laser Doppler flowmeter to 
cold flow model studies. The experimental results to 
date are summarized, and some planned future experi- 
ments are discussed briefly. 


in. BASIC FEASIBILITY EXPERIMENTS 


Since the initial feasibility experiments on the 
laser Doppler flowmeter have been discussed in detail 
in Reference 2, only a brief summary will be presented 
here. The experimental arrangement used to detect 
the Doppler shift in laser light scattered from a flow- 
ing contaminated gas is shown in Figure 2, It is nec- 
essary to add an optical scattering contaminant ( such 
as smoke) to the flowing gas to increase the intensity 
of the scattered light up to a usable level. 



FIGURE 2. DETAILS OF THE EXPERIMENTAL 
"ETUP FOR BASIC FEASIBILITY DEMONSTRATION 

The coherent, monochromatic light beam emerg- 
ing from a He-Negas laser passes through an aperture 
Aj to remove the blue haze which surrounds the red 
laser beam. The beam is then split into two parts by 
the beam splitter Sj, The deflected beam is reflected 
off a flat mirror M lf after which it is attenuated by 
passage through a set of neutral density filters. It 
then passes through u n aperture A3 and a focusing lens 
L3, and is finally deflected into a photomultiplier tube 
by the beam splitter Sj. The rest of tlv original laser 
beam passes through Sj, and is focused to a point by 


the lens Lj. A glass chamber filled with air contam- 
inated with smoke is placed at the focal point of L\. 
The lens L2 collects the light scattered from the smoke 
chamber at an angle 6 with the forward direction of the 
unscattered beam. The aperture Aj limits the total 
solid angle over which the scattered light is inter- 
cepted by Lj. The scattered light collected by Lj 
passes through the beam splitter S2 and enters the 
photomultiplier. The lenses L2 and L3 are positioned 
in such a way tbit the two beams are focused at the 
same spot on the photocathode, and the optical system 
is aligned so that the axes of the two converging beams 
entering the photomultiplier are coincident. 

When the smoke -filled air in the chamber was set 
into motion by means of a small air compressor, a 
high frequency heterodyne signal appeared on the os- 
cilloscope screen. A typical waveform is shown in 
Figure 3, The frequency of the signal agreed approx- 
imately with that to be expected for the Doppler shift 
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Sweep Speed: ZOO p sec/cm 
Sensitivity: ZO mv/cm 


FIGURE 3. DOPPLER WAVEFORM FROM 
MOVING SMOKE AT LOW VELOCITY 

from particles moving with the estimated smoke ve- 
locity, The signal disappeared when the referenci 
beam (reflected from minor M t ) was blocked off 
proving that a true heterodyne signal was being ob- 
served. The signal also disappeared when the smoke 
was at rest. When the flow speed through the chamber 
was varied, the heterodyne signal frequency trackec 
the flow speed linearly. These results verified th 
basic feasibility of a laser Doppler flowmeter for 
measurement of contaminated gas flow velocities. 
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IV. QUANTITATIVE FLOW MEASUREMENTS 
AT LOW VELOCITIES 


It is anticipated that in cold flow model studies, 
velocities up to Mach 2 will be encountered. It is 
therefore necessary to verify th at the laser Doppler 
flowmeter can be used up to hypersonic velocities be- 
fore attempting to apply it directly to a cold flow model. 
Following the initial feasibility demonstration de- 
scribed in the previous section, a program was begun 
to make careful quantitative flow measurements from 
low velocities up to hypersonic velocities. 

The first quantitative flow measurements were 
made with the apparatus shown in Figure 4. The laser 
beam is focused at the center of a 5 mm inside diam- 
eter glass tube through which air is drawn by an ex- 
haust blower. Smoke is added to the air upstream 
from the optics to provide the necessary light scatter- 
ing contaminant. The laser beam emerging from the 



FIGURE 4. EXPERIMENTAL ARRANGEMENT FOR 
LOW-VELOCITY MEASUREMENTS 

glrss tube is intercepted by a lens, attenuated, and 
focused on the photocathode of an RCA-C70042C pho- 
tomultiplier tube. Light scattered at an angle 6 by the 


moving smoke is collected by a lens and also focused 
on the photocathode. The heterodyne signal from the 
photomultiplier is fed into a wideband amplifier, the 
bandpass of which is cut off below 1 kHz by an ex- 
ternal high-pass filter. The Doppler shift frequency 
f D is read out on an analog frequency meter (General 
Radio Type 1142-A). 

The Doppler frequency shift f Q is given in the 
general case by the expression [3] 

fp = ( 2ir) -1 [IT -iT] • v , (2) 

where k and k Q are the wave vectors of the scattered 
and incident beams, respectively, and V is the velocity 
of the scattering particles. Application of equation (2) 
to the scattering geometry shown in Figure 4 gives 



o 


where X = 6328 A is the vacuum wavelength of the in- 
cident laser beam. In deriving equation ( 3) use has 
been made of the fact that k « k , and the directional 
spreads in k and k due to focusing have been ne- 
glected. S 0 

As an independent check on the laser flowmeter 
results, the flow velocity was measured directly by 
timing the passage of the leading edge of the injected 
smoke stream between two points 70 cm apart, one 
upstream and one downstream from the optics. It was 
not possible to measure the flow with high accuracy by 
this method, but the flow speed (~ i cm/sec to ~ 100 
cm/sec) was too small to allow the use of a Pitot tube. 
The frequencies f D measured with the laser Doppler 
flowmeter are plotted as a function of the directly 
measured velocities in Figure 5. A theoretical curve 
calculated from equation (3) is also shown. The re- 
sults are quite satisfactory , especially considering the 
inaccuracies inherent in the direct velocily measure- 
ments. Above a velocity of 100 cm/sec the direct ve- 
locity measurements became too inaccurate to be use- 
ful. The Doppler signal could be tracked, however, 
allthewayup to the maximum velocily v ~ 14. 4 m/sec 
(corresponding to fp ~ MHz which the blower system 
could produce. At heterodyne frequencies higher than 
the 1.5 MHz upper limitof the analog frequency meter, 
the heterodyne signal was monitored on an oscillo- 
scope. 

Auxiliary experiments have been performed which 
indicate that the efficiency of optical heterodyning goes 
through a steep maximum in the immediate neighbor- 
hood of the focal point of the incident laser beam, A 
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FIGURE 5. COMPARISON OF LASER FLOW- 
METER RESULTS WITH THEORY 

similar observation has been made by Kroeger [4] 
using laser light scattered from moving diffuse sur- 
faces. As a consequence of this effect, the Doppler 
shift, which is observed with the laser Doppler flow- 
meter, is due to the flow velocity only in the immedi- 
ate neighborhood of the focal point of the incident 
beam — that is , the laser flowmeter gives a local 
velocity measurement. By moving the focal point 
across the diameter of large flow tubes, velocity pro- 
files have been measured. The flow fields which have 
been used thus far have not been sufficiently laminar, 
however, to allow a good comparison between the 
measured velocity profiles and those predicted by the 
theory of laminar flow. 

Following the low velocity measurements made 
with the arrangement shown in Figure 4, a similar 
optical system was setup on a small wind tunnel in the 
Brown Research Laboratories designed to produce a 
maximum flow velocity of about 50 m/sec. The ex- 
perimental setup is shown in Figure 6. A continuous 
air flow is maintained through the system by an ex- 
haust blower located on the roof of the building. Air 
is drawn into the intake through a set of three porous 
air filters located at the base of the tunnel. Smoke is 
generated in the base region by burning commerical 
smoke candles. The smoke mixes with the incoming 
air and passes through a series of contractions into 
the test section , which has a i" x 3" rectangular cross 
section, A Pitot tube is provided to supply an inde- 
pendent measurement of the flow velocity at the center 
of the test section. Flat optical entrance and exit 


windows for the laser and scattered beams are mounted 
flush with the inside walls of the test sectioh to avoid 
perturbing the flow pattern. Above the test section, 
the wind tunnel cross section gradually expands to a 
10" x 10" square shape. The tunnel is then joined to 
a 10" x 10" vent which connects the tunnel to the ex- 
haust blower on the roof. A bleeder vent is connected 
in parallel with the wind tunnel. By adjusting a control 



FIGURE 6. LOW-VELOCITY WIND TUNNEL 

damper in the bleeder vent, the flow velocity in the 
wind tunnel can be varied over a wide range. 

Velocity measurements on the wind tunnel have 
just begun. The Doppler heterodyne signal has been 
detected from low velocities up to the highest velocity 
(~ 50 m/sec) which the tunnel can produce, but care- 
ful quantitative data have not yet been taken. 


V. PROPOSED FUTURE EXPERIMENTS 


After the completion of the low -velocity wind tun- 
nel work described above , it will be necessary to ex- 
tend the measurements to velocities in the sonic and 
hypersonic range. Preparations are now underway 
to make use of a high speed wind tunnel at Marshall 
Space Flight Center for this purpose. The tunnel in 
question covers the Mach number ranges from 0. 3 to 
0.85 and i.5 to 4.5. The major problem at high 
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velocities'is expected to be the generation and injection 
of sufficient quantities of smoke into the flow stream. 

One possible remedy for the smoke problem is the 
use of a gas laser with considerably higher power than 
that available from the He-Ne laser used in the low 
velocity experiments. As the laser power output in- 
creases, the density of smoke required to maintain a 
given scattered light intensity decreases proportion- 
ately. The possibility of using the recently developed 
Argon gas laser (with an output power of about 330 
mW in a single blue line) is also being investigated. 

Following the hypersonic velocity measurements, 
attention will be turned to the development of a multi- 
dimensional velocity measurement system. Two or- 
thogonal components of local velocity vectors can be 
measured by a trivial extension of the system used to 
measure one component. Although the design of a sys- 
tem to measure three noncoplanar components will be 
somewhat more complex, it does not appear to be im- 
practical. 

The laser Doppler flowmeter should eventually 
find application in a wide variety of gas flow measure- 
ments. Its practical ramifications extend far beyond 
cold -flow model studies, although it was for this pur- 
pose that the present experimental program was begun. 


REFERENCES 


i. Yeh, Y. ; and Cummins, H. Z.: Localized Fluid 
FlowMeasurements withanHe-Ne Laser Spectro- 
meter. App. Phys. Lett., 4, 176, 1964. 


2. Foreman, J. W. , Jr. j George, E. W. ; and Lewis, 
R. D. : Feasibility Study of a Laser Flowmeter for 
Local Velocity Measurements in Gas Flow Fields. 
Tech. NoteR-149, Brown Eng. Co., Inc., Huntsville, 
Ala. , May 1965. (The preliminary results were 
to be published in the August 1 issue of the J. Appl. 
Phys.) . 


3. Lee, K. S. H. : On the Doppler Effect in a Medium. 
Tech. Rep. No. 29, Antenna Lab. , Calif. Inst, of 
Technology, Dec. 1963. 


4. Kroeger, R. D. : Motion Sensing by Optical Heter- 
odyne Doppler Detection from Diffuse Surfaces. 
Proc. IEEE, 53, 211, 1965. 


71 



MEASUREMENT OF GAS TEMPERATURE AND THE RADIATION 
COMPENSATING THERMOCOUPLE 


By 



Glenn E, Daniels 


SUMMARY 


as 


} 


To achieve a uigh degree of accuracy in gas tem- 
perature measurement, a method is devised whereby 
several thermocouples may be wired together to make 
the Radiation Compensating Thermocouple, This 
thermocouple is designed so that the errors of meas- 
urement from the radiation environment will cancel 
out; results of tests, equations, and information are 
presented to aid in its fabrication, , , ^ 


where: T 


I. INTRODUCTION 



= gas temperature 
= surrounding temperature 
= temperature of the sensor 


T 


s 


To evaluate the tests and operation of the Radi- 
ation Compensating Thermocouple described later in 
this report, the first section has been used for a brief 
discussion of the theory of gas temperature measure- 
ment. The available literature covering the theoreti- 
cal problems of heat transfer around a cylinder and the 
derivation of a theoretical expression for the error of 
a thermometer used in measuring a gas temperature 
under the various situations which may occur is 
scarce. Because of the complexity of the theory, the 
actual design of the Radiation Compensating Thermo- 
couple has been based on empirical data obtained by 
laboratory and field measurements. 

The accurate measurement of gas temperatures 
presents many problems because of the large errors 
usually present. An object immersed in a gas reaches 
thermal equilibrium at a temperature which is fre- 
quently different from that of the gas. 

The problem of precise temperature measure- 
ment can best be described as follows: 

The temperatuie which a temperature sensor 
(thermometer) indicates is not usually the tempexa- 
ture of a gas, but a temperature which is a result of 
heat transfer from the gas and its surroundings to the 
sensor by radiation and convection -conduction balanced 
against the rate of heat transfer from the sensor by 
radiation and convection -conduction to the gas and 
surroundings and by conduction from the sensor to its 
support (Fig. i). 


H = heat transfer by convection and 
cv 

conduction between the gas and 
sensor 

H^ s = heat transfer by radiation between 
? the sensor and the surroundings 

**CD = heat transfer by conduction 

between the sensor and support 

Hrq = heat transfer by radiation from 
gas molecul to sensor 


FIGURE i. HEAT TRANSFERS BETWEEN AN 
OBJECT AND GAS 


If the temperature of the gas changes rapidly, the 
response time of the temperature sensor (i. e. , the 
time required for the individual heat transfer rates to 
again come into equilibrium) is another factor which 
influences the correctness of the temperature indi- 
cated by the sensor with respect to the true gas tem- 
perature. 

The knowledge that thermometers do not usually 
measure the correct temperature of a gas is not new; 
William Charles Wells [i] in his studies of the forma- 
tion of dew at the beginning of the 19th Century ( 1800 
to 1815) realized that the air temperatures measured 
by his thermometers were not always correct. The 
Annual Report of the Chief Signal Officer [2] for the 
year 1887 includes a detailed discussion of the errors 
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of thermometers used to measure air temperatures 
and of the errors of air temperature measurements 
when weather shelters are used. 

Since these early reports were made, a large 
number of studies and tests on the problems of gas 
temperature measurement have been reported in the 
literature; for example, Circular 513 [3] of the Na- 
tional Bureau of Standards lists approximately 400 
references to papers on the measurement of gas tem- 
peratures published during the period between 1929 
and 1951. Since 1951, many additional reports have 
been published on the subject. 

In many published articles, the factors which 
cause a sensor to indicate a temperature different 
from the gas temperature have been described as 
errors of the sensor. Actually, the temperature indi- 
cated by the sensor is a very accurate measurement 
of the equilibrium temperature the sensor has reached 
while exposed to the gas and to its surroundings, if 
properly calibrated. 

The amount by which the temperature of the 
sensor differs from the actual temperature of the gas 
is a result of the various heat transfers reaching a 
balance different from the balance of the heat trans- 
fers which would result if the temperature of all the 
surroundings were the same as the gas temperatu/e. 
This amount will be referred to as an error in this 
paper (not an error of the sensor, but a result of the 
sensor temperature being different from the gas tem- 
perature) . 

The various factors which enter into the meas- 
urement of a temperature affect a temperature sensor 
in a gas in the following ways: 

1. Radiation heat transfer (H^g) between the 
sensor and the surrounding objects occurs when the 
surroundings have a different temperature from the 
gas; i. e. , the sun, as a source hotter than the sensor, 
will cause the unshielded thermometer to indicate too 
high a temperature, while a clear sky, as a source 
colder than the sensor, will result in temperature in- 
dications which are too low. The addition of thermal 
radiation shields reduces, but does not eliminate, the 
radiation error because the shield (a surrounding 
object) becomes a source for heat transfer by radia- 
tion. 

Reduction of the amount of heat transfer by radia- 
tion (i. e, , a cloudy day) will reduce the error of 
measurement. Radiation heat transfer can occur in a 
vacuum. 


% 

2. Radiation heat transfer between the sensor 
and the surrounding gases (H^q) is an error neg- 
ligible in meteorological measurements because air is 
transparent to radiation transfer; but luminous parti- 
cles in gases in furnaces, f rom certain fuels, cause 
errors the same way that solid objects cause radiation 
errors. 

3. Convection-conduction heat transfer between 
the sensor and the gas (H^y) is the main process 
through which the temperature sensor approaches the 
gas temperature. If this were the only method of heat 
transfer, the sensor would always reach the gas tem- 
perature very closely if its response time vere fast 
enough. This heat transfer will vary with the speed 
of motion of the gas across the sensor, being lowest 
with calm winds and increasing as a function of the 
wind speed. The greater this heat transfer the faster 
and closer the temperature sensor will approach the 
desired gas temperature. 

4. Conduction heat transfer from the sensor to 
the sensor supports (Hq d ) is a loss of heat (or 
gain) that could cause large errors. This is the flow 
of heat from a thermocouple along the wire leads, or 
the flow of heat from a thermometer bulb to the glass 
stem. This heat transfer increases as the thermal 
conductivity increases, but decreases with decrease 
of the cross section of the supports. Such heat trans- 
fer needs to be small for accurate temperature meas- 
urement. 

5. The impact of gas molecules on a sensor can 
cause an error by heating the sensor. For low gas 
velocities (normal wind speeds), these errors are 
negligible, but at high velocities, such as space ve- 
hicle velocities, these errors can be important. 

6. The response time of a sensor can be a factor 
which introduces an error in measurement when a 
sudden change in gas temperature occurs. The faster 
the sensor response, the smaller the magnitude of 
this error and the shorter the time in which it occurs. 

There are several ways in which the errors can 
be reduced: 

1. When the sensor size is reduced, the heat 
transfer per unit area from radiation will be the same. 
The response time will also be shortened. The heat 
transfer per unit area from conduction and convection 
to the gas will increase. Therefore, the smaller sen- 
sor will more closely approach the true gas tempera- 
ture; i.e. , it will have a smaller error. Margaret 
Fishenden and Owen A. Saunders [4] have given a 
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comprehensive discussion of heattrans'fer on cylinders 
of various sizes. This reference should be consulted 
if more information is desired. 

W. J. Waggener [ 5] , at a scientific meeting in 
Germany in 1898, proposed, probably for the first 
time, the use of a group of thermocouples of various 
wire sizes to determine the true temperature of a 
Bunsen burner flame ( a hot gas surrounded by colder 
objects). His method was to place several thermo- 
couples of different sizes in the flame; the various 
temperatures obtained were plotted versus the size of 
the thermocouple wires. By extrapolating the tem- 
perature curve obtained, to zero wire size, the true 
gas temperature is determined, because a thermo- 
couple of zero mass and size will theoretically indicate 
the correct gas temperature. This is illustrated in 
Figure 2. 



FIGURE 2. USE OF THERMOCOUPLES OF 
VARIOUS SIZES TO FIND THE GAS 
TEMPERATURE, T 

o 


Henry Kreisinger and J. K. Barkley [ 6] conducted 
a series of tests, using thermocouples of various 
sizes to measure temperatures of gases in boilers. 

They found that when measuring temperatures of hot 
gases surrounded by cooler surfaces the temperatures 
measured by small thermocouples were higher than 
those of large thermocouples. They also found that 
the higher the gas temperature in relation to the tem- 


perature of the boiler walls, the greater the error of 
measurement of the true gas temperature for a sensdr 
of a specific size. When using a thermocouple 0, 020 
cm in diameter, errors exceeding 55° C (100®F) were 
found when the gas temperature was about 1090° C 
(2000®F) . 

Data collected at Fort Huachuca, using copper 
and constantan butt-joined thermocouples of various 
sizes in measurement of air temperatures, have 
shown that, when the thermocouples are used without 
shields or special reflective coatings, errors can 
reach 1°C for AWG No. 35 gage wire (0. 014 cm di- 
ameter) during clear days and one-half this value on 
clear nights. If the thermocouple is reduced in size 
sufficient to provide accuracies of air temperature 
measurements of ±0. 06° C (0. 1°F) , as stated in many 
requirements for accuracy of air temperatures, a 
thermocouple no larger than 0. 0008 cm diameter 
(AWG wire size about No. 58) will be required. Such 
a thermocouple would be difficult to fabricate, install, 
and maintain and could not be expected to remain in- 
tact for more than a short period of time. 

2. If the surface emissivity of the sensor is re- 
duced (reflectivity increased) , the absorption of ra- 
diation will be lower (heat transfer from radiation is 
less). Therefore, the sensor will read closer to the 
temperature of the gas. One difficulty with using high 
reflectivity coatings is that many of them have high 
values of reflectivity in the visible region of the spec- 
trum (appear white) , yet in the infrared region of the 
spectrum they have a low reflectivity, and thus act as 
a black body [ 7] . 

Considerable study of the problem of coatings for 
temperature sensors for radiosondes has been made 
by the University of Michigan and published in an 
annual progress report [ 7] . One of the best coatings 
they found, with a high reflectivity in both the visible 
and Infrared region of the spectrum, is aluminum. It 
has a reflectivity of 94 percent. Such a coating on a 
thermocouple made of AWG No. 35 gage wire would 
result in the error of measurement, because of radia- 
tion, of about 0. 17° C (0. 3®F) when used in the sun 
without shields. This error is much smaller than 
normal measurements, and continuous measurements 
of air temperatures for micrometeorological use with 
this small error would be of great value. The dis- 
advantage of such a thermocouple is that a very small 
thermocouple is required, and even small particles oi 
dust accumulating on the thermocouple will greatly 
increase the measurement errors. 

3. Since increasing the speed of the gas over the 
sensor increases the heat transfer from conduction 
and convection in the gas, aspiration of sensors has 
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long been a favorite method of reducing the error of 
gas temperature measurement. The "high velocity" 
thermocouple is most useful in measurement of gases 
in furnaces where errors of 200°C are common on the 
large thermocouples, where an appropriate design of 
the sensor and shields to accelerate the flow of gas is 
used. The "high velocity" thermocouple when thus 
used can be placed parallel to the flow of the gases. 
In meteorological measurements of air temperatures, 
especially near the surface of the ground, aspiration 
will obtain air samples, not from a point or horizontal 
plane outside the shield but from a sphere of air, 
with its radius depending on the rate of aspiration. 
Wind, wi*h its changing direction and velocity, will 
distort the sphere so that the center of the sphere will 
not coincide with the location (height above the ground) 
where the measurement is required. J. Bartel [ 8] 
and Alf Nyberg [ 9] have discussed this problem in 
detail. 

Many methods have been reported in the literature 
to eliminate, reduce, or correct the error by which 
sensors differ from true gas temperatures. E. Marston 
Moffatt, [10] Andrew I. Dahl, [11, 15] W. L. Bolies, 
[12] W. L, Servinghouse, [13] Marvin D, Scadron, 
Isidore Warshawsky, and Clarence C. Gettelman, 

[ 14] have all discussed these gas temperature errors 
in detail and give methods for estimating the errors. 
The National Advisory Committee for Aeronatuics 
(now Lewis Research Center, National Aeronautics 
and Space Administration) conducted studies on this 
problem, to use in research on jet engines, and pre- 
sented the results in Technical Reports [16, 17, 18] . 
In NACA Technical Note 3766 [ 17] the following equa- 
tion was given for computing the error of bare-wire 
thermocouples: 


If we examine equation (1), we see that many of 
the parameters in the equation would be the same for 
two thermocouples of different size but of identical 
construction (type, same materials, etc. ) , if they were 
exposed close together. Therefore, we can rewrite 
equation (1) as follows for forced convection: 


Radiation correction = function 
where 



(la) 


d = diameter of thermocouple wire 
v = velocity of gas across the sensor. 

Equation (la) shows that the smaller the sensor 
the smaller the error; likewise, the higher the aspi- 
ration rate the smaller the error. For a short enough 
period of time the velocity of the gas could be assumed 
constant and then the radiation correction is a func- 
tion only of wire size. 

If the gas is at a temperature, T 0 , and the sur- 
rounding objects are at a different temperature, T s , 
(represented as the sum of all the objects at tempera- 
tures T Sj , Tg 2 , T S3 T Sn ) lower (colder) than 

T 0 , we can illustrate the various heat transfers and 
the flow of heat as shown in Figure 1. 


The condition for equilibrium of the sensor in 
Figure 1 can be shown by 


H cv + h rg + h rs + h cd 0 


( 2 ) 


where 


Radiation Correction = 


K» rad 



-o- 18 



where 


= heat transfer by convection-conduction between 
the gas and the sensor. 

- heat transfer between the gas molecules and 
the sensor by radiation 


K* . * 27e 
rad w 

T = indicated thermocouple junction temperature 
= 1000° R, reference temperature 
= temperature of duct (surroundings) 

M 0 - velocity of gas (in terms of Mach number) 

r 

e = emittance of wire 
w 

d = diameter of wire 


H rs = heat transfer by radiation between the sensor 
and the surroundings. 

H^ d = heat transfer by conduction between the sensor 
support and the sensor. 

Some terms are negative or positive depending on 
whether surrounding objects are warmer or colder 
than the gas, and whether the sensor is warmer or 
colder than the gas or the supports. If the gas is 
transparent to radiation (as is air), H^q is zero, 
then we can write equation (2) as 

b 0 <T 0 - T t ) A, + «se 8 (T| - T{) A,X + <KA 2 /a)(T t - T a ) =0, (3) 
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where 

h = coefficient of heat transfer by convection- 
conduction 

T q = gas temperature 

T^, = sensor temperature 

T g = temperature of surrounding objects 

A t = surface area of sensor 

A 2 * cross-sectional area of sensor 

a = Steffan- Boltzmann constant 

e = surface emissivity of sensor 

e = emissivity of surrounding objects 
s 

X = fraction of total solid angle over which the junction 
sees the surrounding objects 

K = coefficient of thermal conductivity of sensor and 
support 

a = depth of immersion of sensor and support in the 
gas. 

Tjie term aee s (T| - T T ) AjX will nearly always 
consist. of a sum of several terms, one term for each 
variety of surrounding source of radiation, of tem- 
perature, or of emissivity. The term for conduction 
(KA/a) (T>j> - T s ) likewise may consist of several 
terms, one for each type of material used in the sen- 
sor or support. 

Equations (1) and (2) show that, when the net 
temperature of the surroundings differs from the gas 
temperature, a sensor will indicate not the tempera- 
ture of the gas, but some temperature between the gas 
temperature and the surrounding temperatures. 

Therefore, theoretically a single temperature 
sensor never reads the true temperature of a gas ex- 
cept in a special case where no heat transfer occurs. 
This fact is important because we do not have a direct 
method of precisely measuring the correct tempera- 
ture of a gas and therefore have no reference to use 
to determine directly the error of any gas thermometer 
from the true gas temperature. Too often experiments 
state that their thermometers have no error when they 
actually have no way of knowing the errors of their 
measurements. But it is possible, by proper selection 
of several sensors, identical except for size and 
properly connected such that the error on each sen- 
sor when summed up with proper sign will cancel out, 
that the gas temperature can be measured directly 


with a negligible error. This is the principle of the 
Radiation Compensating Thermocouple. Its design is 
such that it can be connected directly in place of any 
single thermocouple of like materials, but will give 
measurements directly of gas temperatures with much 
smaller errors. Any shields or aspiration required 
on the single thermocouple can be eliminated. 

II. THE RADIATION COMPENSATING 
THERMOCOUPLE 

If we look atFigure 3, which is Figure 2 redrawn 
for a condition where the gas temperature (air tem- 
perature) is lower than the surrounding source (sun, 



FIGURE 3. THEORY OF RADIATION COMPENSAT- 
ING THERMOCOUPLE - STATIC CONDITION 

on a clear day plus the sky) , we have a curve of tem- 
perature where the temperature increases with sensor 
size (thermocouple wire size). The extrapolation of 
the curve to zero sensor size gives the true air tem- 
perature (T 0 ). In this case we are considering one 
individual static case of an unchanging air tempera- 
ture and surrounding temperatures such as that which 
occurs for a short time on a clear day. For this con- 
dition, we can select three sensors of three sizes 
(small, medium, and large). The temperatures 
measured by these three sensors can be represented 
by T s , Tjyj, and T L , respectively, with errors AT S , 
AT,^, and ATl from the true gas temperature. 
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Therefore, 

T = T + AT 
so s 


T = T + AT, . 
Mo M 


T. =T + AT. . 
L o L 


The three thermocouples may be wired together 
such that thermocouple Tl is opposing the other 
thermocouples electrically. Then these three ther- 
mocouples will read some temperature T x . and we 
can write 

T * ■ ,T „ + at m> ■ < T ,, ♦ AT L» + < T „ + AT „> • < 5 “> 

or, simplifying, 


T = T + AT,. - (AT t - AT ). 
x o M ' L s' 


If we select these three thermocouples such that 
the errors satisfy the equation 

aT M * AT L - 4T S ' <«> 



o 


FIGURE 4. WIRING OF RADIATION COMPENSAT- 
ING THERMOCOUPLE 


then from equation ( 5b) we have 

T = T . (7) 

x o ' ' 

,ther words, by proper selection of the three 
thermocouples, the system can be made to read the 
true gas temperature directly. The wiring of such a 
system is shown in Figure 4. 

The previous discussion is based on a single static 
condition of radiation error. If we look at equation 
(la) we see that the radiation correction (radiation 
error) is a function of the square root of the diameter 
of the thermocouple and therefore a change of radia- 
tion from the surroundings, change of gas temperature, 
or change of wind speed to new constant values (for 
new static conditions) will shift the curve in Figure 4 
to a new position with the same mathematical shape. 
The individual errors of each single thermocouple 
may be smaller or larger, but because the intersection 
of the new curves by the same wire sizes will be by 
parallel lines, the ratio 


M 8 

AT. - AT (8) 

Xi s 

will remain a constant and equal to "one" for the ideal 
thermocouple design. 


The theoretical description of the Radiation Com- 
pensating Thermocouple has been simple, but the 
actual design of a working unit, i. e. , finding the cor- 
rect wire oizes, required careful and precise tests 
to evaluate all phases of construction. 

The Radiation Compensating Thermocouple, when 
evaluated in terms of the errors discussed earlier, 
gives the following: 

1. The error, from heat transfer by radiation 
between the sensor and the surroundings, and the heat 
conduction-convection between the gas and sensor not 
being in biuance, can be reduced to any accuracy de- 
sired depending on the precision of design and fabri- 
cation. 

2. The error from heat conduction between the 
sensor and support can be considered in the design. 
This error is small enough to be ignored, with the 
present design. 

3. The Radiation Compensating Thermocouple 
does not compensate for response time error; In fact, 
if the thermocouples selected vary too greatly in size, 
the compensating balance will change for short periods 
of time with rapid temperature changes. It was for 
this reason that the original design using three therm- 
ocouples was changed to one using five thermocouples. 
With five thermocouples , we can rewrite equation ( 6) as 
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(84) 


TABLE 1. 


But, If the two large and two small thermocouple© are 
the same size, 

AT., = 2 (AT - AT ) (6b) 

MLs 

or 

? iT « ,Ai i- iT 5 - <fc) 

This permits fabrication with thermocouples more 
nearly the same size and reduces the effect of re- 
sponse time error. 

lit RESULTS OF TEST IN THE METEORO- 
LOGICAL REGIME 

The testa which led to the discovery and were 
used in development of the Radiation Compete ating 
Thermocouple were conducted in the atmosphere in 
the temperature regime used for meteorological meas- 
urements. The equipment available* for recording nod 
for tests strongly influenced the physical design of the 
prototype models. Copper-cooatantan type thermo- 
couples were used, mounted on a tower with the ther- 
mocouples facing south. This required a design of the 
Radiation Compensating Thermocouple with the five 
Junctions and the wire leads connecting the Junctions 
forming a flat surface parallel to the ground (Fig. £ 
and 6). The five junctions were along the south-facing 
edge* of the plane and theleads horizontal. In this way, 
all the thermocouple Junctions had similar exposure to 
the sun (the largest contribution of point type radia- 
tion) . If the thermocouple is rotated around Urn north- 
south axis, there will be little change 1ft the compen- 
sating ability except for convection between individual 
thermocouples unequally affecting other thermocouple* 
in the thermocouple system. This convection would 
exist In its worst condition when the plane of the 
thermocouple system Is vertical since the lowest 
thermoooupu receives no convection from the other 
thermocouples. Yet aU of tha other four thermocouples 
would be affected by different amounts of convection. 

All teste conducted were controlled by five indi- 
vidual thermocouples to give a reference true air tem- 
perature. Except for size, these thermocouples were 
as identical, in construction as possible (Table 1). 

To obtain the. largest possible errors on each 
single thermocouple, each thermocouple element in 
the Radiation Compensating Thermocouple was given 
a thin coating of flat black paint (Kodak black brush- 
ing L&oqueryrajs used) to make the absorption of radia- 
tion as great a« pqeaible. 


DIMENSIONS OF SINGLE THERMOCOUPLES USED 
IN TESTS 


I am of Taamoomw 

AUG 

9ir« files 

tf&f« 5 (m t*r* 

CJschsa) 

VUtfc ajf Thtrncemipl* f*C*«u 

(luefru) 

35 

».<W 

0.025 

21 

0.011 

o.oio 

34 

0,020 

0.100 

20 

D,Q« 

0.1W 

11 

o.on 

0.110 


* Tbexwi cot^let tlMfftvn AXtAtdeJ cHaaa *ilu*a by 


than 0.003 inch rejcccAd. 



FIGURE 6. INSTALLATION OF RADIATION 
COMPENSATING THERMOCOUPLE 


In the early tests at the Meteorological Depart- 
ment, Fort Hu ac hue ft, Arizona, the only recording 
potentiometer* available for data collection were 
multiple print recorders with a printing cycle of 15 
seconds between aach print. Because of this, one 
minute and 15 seconds was required to obtain tem- 
perature readings from the five single thermocouple* 
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FIGURE 6, MODEL OF RADIATION COMPENSATING 
THERMOCOUPLE 

and one Radiation Compensating Thermocouple. Much 
of the data collected during the day time when the sun 
was out was unsatisfactory because of turbulence 
which frequently caused air temperature changes of as 
much as 10° F in afew seconds during each test cycle. 
The early morning data did provide short periods of 
data that could be used in design of the early test 
models. Some useable data were also obtained by 
watching the recorder and determining the tempera- 
ture difference between each reading. Later a syn- 
croprint recorder (with printing speed usually less 
than one second per point) was available. The data 
obtained in the morning hours with chis recorder were 
frequently excellent, and by studying the records, 
periods of high stability and calm winds could botound 
which had repeated cycles of data where each thermo- 
couple would read nearly the same temperature dur- 
ing successive printouts. 

Data collected on one of the tests with the syncro- 
print recorder was selected as the best data for use 
in designing models of the Radiation Compensating 
Thermocouple for future tests. The data given in 
Table 2 are values averaged over several cycles of 
the recorder and represent data where the radiation 
error on each single thermocouple was the largest 
observed. The use of this data for design of the 
Radiation Compensating Thermocouple for regular use 
in measurement of air temperatures would mean that 
most of the time the errors on each single thermo- 
couple would be less than that used for design? there- 
fore, theRadiation Compensating Thermocouple would 
usually be performing better than the design data in- 
dicated. 

Using the data in Table 2 and assuming a rela- 
tionship of 

T = a + b D + c D 2 (9) 

where 

T = air temperature 


D = (diameter of wire) , 

constants a, b, c, were determined by least squares 
solution and the following equation was obtained: 

T = 69. 112 + 279. 832D - 1033. 06D 2 (10) 

If we substitute a wire diameter of zero into the 
equation, then we find a value for T Q (zero wire size) 
of 69. ll°F;i.e. , the air temperature was 69. 11°F at 
the time. 

TABLE 2. 

AVERAGE VALUES OF TEMPERATURE RECORDED 
AT FORT HUACHUCA, ARIZONA, 09:00A. M. APRIL 
8, 1960, BY VARIOUS SIZE THERMOCOUPLES (ALL 
THERMOCOUPLES PAINTED BLACK AN”D EXPOSED 
WITHOUT SHIELDS) 


Thermocouple 1 

Average Temperature 
CF) 

Wire Size 
(AWG Size) 

Diameter** 

16 

0.053 

81.06 

20 

0.0335 

77.26 

24 

0.022 

74.79 

28 

0.014 

72.91 

35 

0.006 

70.70 

Model 42 1 

Radiation Compensating Thermocouple 

68.91 


** Actual wire diameter measured after painting. 


Now, let 
AT = T - T 

o 

v/here AT = temperature difference from T 0 . Sub- 
stituting into equation (10) , we obtain 

AT = 279. 832D - 1033. 06D 2 for AT in °F, and D 
in inches (Ha) 

or 

AT = 6. 1205D - 0. 88958D 2 for AT in °C, and D 
in millimeters. (lib) 

Equations (11a) and (lib) are empirical equations 
which canb' used to find the error of any small black 
temperature sensor in the sun. See Table 3 for errors 
of thermocouples computed from equation (11a) and 
shown in Figures 7 and 8. If the error ; s desired for 
a sensor which has an emittance different from a black 
body, the error can be found by correcting the error 
obtained from equation (11a) or (lib) by the value of 
the emittance. 
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TABLE 3. 

TEMPERATURE ERRORS OF SINGLE THERMO- 
COUPLES OBTAINED FROM DATA IN TABLE 2 
EQUATION (iia) COMPARED TO ACTUAL 
MEASURED ERROR IF T 0 = 69. 11 

COMPUTED FROM 

/&T( 0 F) = 279.83 D (INCHES) 1033.06 D 2 (INCHES 2 ), 

OR 

AT(°C) = 6.1205 D (mm) - 0.88958 D 2 (mm 2 ) 


Wire Size (D) 

AT 

AW6 

(INCH) 

(mm) 

°F 

°c 

16 

0.0508 

1.2903 

11.54 

6.41 

18 

0.0403 

1.0236 

9.60 

5.33 

20 

0.032C 

0.8128 

7.89 

4.38 

22 

0.0254 

0.6452 

6.43 

3.57 

24 

0.0201 

0.5105 

5.20 

2.89 

26 

0.0159 

0.4039 

4.19 

2 # 33 

28 

0.0126 

0.3200 

3.36 

1.87 

30 

0.0100 

0.2540 

2.70 

1.50 

32 

0.0080 

0.2032 

2.17 

1.21 

34 

0.0063 

0.1600 

1.72 

0.96 

35 

0.0060 

0.1524 

1.64 

0.91 

36 

0.0050 

0.1270 

1.37 

0.76 

50 

0.0010 

0.0254 

0.28 

0.16 


The actual errors of the single thermocouples can 
be compared to the errors computed from equation 
(11a) , as shown in Table 4. 


the measuring thermocouple and the small thermo- 
couple sizes were found such that the measuring 
thermocouple was approximately halfway between the 
large and small thermocouples on the error curves. 


TABLE 4. 

ERRORS OF SINGLE THERMOCOUPLE SENSORS 
(BLACK, BUTT-JOINED) FOR VARIOUS WIRE 
SIZES, EXPOSED IN BRIGHT SUNLIGHT 


| Thermocouple 




Hire Sire 
AWG 

Diameter 

(Inches) 

Computed 

•P 

Actual 

•F 

Computed Minus Actual 
•F 

16 

0.053 

11.93 

11.95 

-0.02 

20 

0.0335 

8.21 

8.15 

+0.06 

24 

0.Q22 

5.65 

5.68 

-0.03 

28 

0.014 

3.72 

3.60 

-0.08 

35 

0.006 

1.60 

1.59 

+0.01 


NOTE: Model 42 of the Radiation Compensating Thermocouple (data la 
Table 2) had an error of 0.20*F from the computed T 0 . Thi* 
la within the measuring accuracy of the equipment being used. 



To construct a Radiation Compensating Thermo- 
couple, equation (11a) or (lib) can be used to deter- 
mine the wire combinations. The method used in the 
tests was to select a specific size thermocouple for 
the largest thermocouple; and then, by trial and error, 


FIGURE 7. ERRORS OF SINGLE THERMOCOUPLE 
SENSORS (BLACK, BUTT- JOINED) FOR 
VARIOUS WIRE SIZES, U. S. CUSTOMARY 
UNITS 
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1.4 


DESIGN II 



Thermocouple 

Wire Size 




(A WG) 

(INCHES) 

AT 

^L*^S 

S 

28 

0.0126 

3.63 

4.53 

L 

20 

0.032 

7.89 


M 

24 

0.0201 

5.20 


L 

20 

0.032 

7.89 

3.70 

S 

26 

0.0159 

4.19 



£ <^- 

^S> 


8.23 


Using equation (8), we obtain 


AT 

M _ 5.20 

AT. - AT 8. 23 
Li o 


0.63 . 


Loo 


Itemperaftre Error (AT) *C 


FIGURE 8. ERRORS OF SINGLE THERMOCOUPLE 
SENSORS (BLACK, BUTT-JOINED) FOR 
VARIOUS WIRE SIZES, METRIC UNITS 

As examples of verification of the amount of error 
in specific designs of the Radiation Compensating 
Thermocouples, the following designs are analyzed 
with equations (11a), (6a), and (8). 


In Design II, the Radiation Compensating Therm- 
ocouple is overcorrectiug the error by about 37 per- 
cent. This means that the Radiation Compensating 
Thermocouple will read a temperature lower than the 
actual gas temperature when a radiation source hotter 
than the surroundings is present; i. e. , the tempera- 
ture indicated by the Radiation Compensating Therm- 
ocouple will be higher when shaded than when in the 
sun. 


DESIGN El 


DESIGN I 


Thermocouple 

Wire Size 

AT 

5.20 

7.89 

6.43 

7.89 

5.20 

- AT S 
2.69 

2.69 

S 

L 

M 

L 

S 

■ . ■ 


e - a. r s > 

5.38 


Using equation (8) , we obtain 
AT. 


M 


_ 6,43 

AT t - AT., 5.38 

Li o 


= 1 . 20 . 


This means the Design I has a residual error of about 
20 percent still in the measurement. In other words , 
it has removed 80 percent of the error. 


Thermocouple 

Wire Size 


^L‘^S 

1.84 

2.50 

S 

L 

M 

L 

o 

(AWG) 

28 

24 

26 

24 

30 

(INCHES) 

0.0126 

0.0201 

0.0159 

0.0201 

0.0100 

A}? 

3.36 

5.20 

4.19 

5.20 
2.70 

Z - £Tg) 

4.34 


From equation (8), we obtain 



AT L ~ AT s 


4. 19 
4. 34 


0.97 . 


Design III will provide a Radiation Compensating 
Thermocouple that will have approximately 3 percent 
error in measurement of the true gas temperature 
indicated being lower than the actual gas temperature. 
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All computations and discussions of the designs 
were based on the worst conditions: a clear, sunny, 
calm day and black thermocouples. To get an idea of 
the Radiation Compensating Thermocouple compared 
to a single thermocouple under similar conditions, let 
us examine Design III. The single thermocouple, 
designated "M," would normally be used in measure- 
ment. If all thermocouples are black, the n M M therm- 
ocouple would have an error of 4. 19° F on a clear day 
in the sun. At the same time, the Radiation Compen- 
sating Thermocouple would have an error of 4. 19 - 
4. 34 * - 0, 15°F; i.e. , it would read 0. 15° F too low. 

This is the '’worst” error to expect. Also, these 
computations were made on the assumption that the 
fabrication would not increase any of the diameters of 
the wires. Actually, about 0. 001 inch should be added 
to each wire size before making computations. 

Design II has provided the best data in tests be- 
cause of the adjustment in diameters of the wires for 
fabrication. 

These tests were carried on over a period of about 
two years. During this time, about 25 models of the 
Radiation Compensating Thermocouple were fabricat- 
ed by hand using ordinary jewelry tools, pliers, small 
files, silver solder and a propane torch. Besides 
comparing the Radiation Compensating Thermocouple 
with single thermocouples, two identical Radiation 
Compensating Thermocouples, designed to give very 
small errors, were fabricated. One was painted 
black, and the other white using a flat white paint that 
was white not only in the visible region of the spec- 
trum but also in the infrared region of the spectrum. 
When these two thermocouples were placed in the sun, 
the air temperatures recorded were almost exactly 
the same from each Radiation Compensating Thermo- 
couple under stable air conditions. This is an excel- 
lent example of the usefulness of the Radiation Com- 
pensating Thermocouple in gas measurements, since 
contamination or changes in the surface color do not 
change the accuracy of the Radiation Compensating 
Thermocouple unless the contamination is not uniform 
and thus changes the balance between the individual 
thermocouples. 

Since the original tests, operational tests have 
been made at Fort Huachuca by Mr. A. V. Carlson of 
the Meteorology Department [ 19] . In these tests 
various test models of the Radiation Compensating 
Thermocouple were compared with a single AWG No. 
24 gage thermocouple. These tests were recorded on 
a two-pen recorder permitting simultaneous record- 
ing of temperature data from the single thermocouple 
and the Radiation Compensating Thermocouple, These 


tests were run with various radiation conditions - 
clear and cloudy days and nights with the thermo- 
couples inside and outside the standard instrument 
shelter. 

The results of these tests were as follows: 

1. The performance of the test models of the 
Radiation Compensating Thermocouple was as pre- 
dicted. (These predictions were based on the data 
given in this report . ) 

2. A properly designed Radiation Compensating 
Thermocouple used without shields or shelter is su- 
perior to a single thermocouple even when the single 
thermocouple is mounted in an instrument shelter. 

3. The single thermocouple read 4 to 5 degrees 
Fahrenheit higher than the Radiation Compensating 
Thermocouple when both were exposed in direct sun- 
shine even in winds up to 15 miles per hour. 

One more remark is needed on the data in Table 
2. These data are only one sample from one morning; 
but this is the best set of data obtained for the worst 
conditions during the two years of tests. All other 
data agreed closely with this data within the errors of 
data measurement, and because of the excellent fit 
of this data by a second-degree curve, these data 
should provide adequate information for preliminaiy 
design of Radiation Compensating Thermocouples 
within the thermocouple wire size range given and 
general meteorological temperature range. The tests 
at Fort Huachuca [19] have further substantiated the 
validity of this data. For temperatures outside the 
meteorological range, tests should be made with five 
or more single thermocouples of various sizes in the 
desired temperature ranges to verify or correct the 
curves. 

As a check for faulty construction, all Radiation 
Compensating Thermocouples fabricated should have 
a laboratory or field test to compare them with the 
extrapolated temperature for wire size zero, using a 
group of single thermocouples. This could be done in 
a laboratory test chamber. 

Patent 3, 049, 012 on this device was issued August 
14, 1962, with rights assigned to the United States 
Government., inventor Glenn E. Daniels. 

IV. EFFECTS OF PRESSURE CHANGES 

The NACA equation, discussed in this report as 
equation ( i) , for computing the radiation error of a 
thermocouple does not include any term for gas pres- 
sure, and is valid only at standard atmospheric pres- 
sure. Actually, the radiation error increases with 
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decreasing atmospheric pressure and becomes very 
large at low pressures. This is because the major 
factor which opposes and reduces the radiation er- 
ror, heat transfer convection-conduction between the 
gas and sensor, becomes smaller with decreasing air 
pressure, because there are fewer molecules for the 
convection and conduction. The values of error with 
sensor size, which are given in Table IV, were ob- 
tained at Ft. Huachuca, Arizona* The average sta- 
tion pressure at Ft* Huachuca is about 0, 8 of an 
atmosphere; therefore, the measured temperature 
errors in the tests at Ft. Huachuca will be a little 
larger than at sea level, but smaller than higher ele- 
vation stations. 

V. CONCLUSIONS 


Current gas temperature measurements have 
large errors. The Radiation Compensating Thermo- 
couple can provide data on the temperature of a gas 
much more representative of the actual gas tempera- 
ture than either a single thermocouple, with or with- 
out shields, or an aspirated thermocouple. Even if 
the Radiation Compensating Thermocouple is crudely 
made, and the error is not completely compensated 
out, the data from this thermocouple are superior to 
data collected by the sensors currently in use. 

Additional tests should be made in the laboratory 
and in the field on the Radiation Compensating Therm- 
ocouple to better understand the theory of heat trans- 
fer and to determine the optimum design for its use 
in various types of gas measurements. 

Tests made at Ft. Huachuca, Arizona, proved the 
feasibility of the Radiation Compensating Thermo- 
couple in gas temperature measurement. 

For immediate use, the information given in this 
report can be used to design Radiation Compensating 
Thermocouples which have small errors in gas tem- 
perature measurements when used without shields. 
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SUMMARY 

.. 29 3 -'A'- 

Statistical analyses of turbulent fluctuations are 
needed to provide inputs for dynamical and structural 
response calculations. Previous experimental inves- 
tigations of area-integrated correlation functions are 
limited because they were based on point measure- 
ments with standard solid probes like hot wires and 
dynamic pressure transducers. In this paper the use 
of the integrating features of optical beams for a "one 
shot" estimate of area integrated correlation functions 
is discussed. ^ * * 

INTRODUCTION 


The statistical analysis of turbulent fluctuations 
is needed to provide inputs like aerodynamic forcing 
functions, turbulent fluxes, and Mach wave sound 
sources for dynamical and/or structural response cal- 
culations. The mathematical description of all these 
inputs is identical insofar as they are expressed by an 
area integral over a space -time correlation function. 
The same integrals are needed in the statistical theory 
of turbulence, if one wants to calculate true three- 
dimensional wave number components instead of their 
one -dimensional approximation based on line integrals 
and symmetry assumptions fl]. 

Previous experimental investigations of area in- 
tegrated correlation functions are limited because they 
were based on point measurements with standard solid 
probes like hot wires and dynamic pressure trans- 
ducers, The inherent numerical integration proce- 
dures require a large number of points. Ihe experi- 
mental and numerical effort in providing these 
pointwise estimates of correlation functions is then 
prohibitive, even if an in cf ^ument could be found that 
combines the necessary temporal and spatial resolu- 
tion with a sufficiently linear and time invariant fre- 
quency response [2]. 

We propose to use the integrating features of op- 
tical beams for a "one shot estimate" of area integrated 
correlation functions. The cross correlation of light 
incensities gives the three -dimensional wave number 
components of light extinction coefficients for any wave 
front propagating normal to the plane of the crossed 


optical beams. The wave-length and the spectroscopic 
resolution of the received light may be adjusted to 
represent a selected thermodynamic properly such as 
density. An application to Mach wave sound emission 
leads to a one shot estimate of Ffows Williams* Mach 
wave sound sources [3], 

I. WAVE NUMBER COMPONENTS OF LIGHT 
EXTINCTION COEFFICIENTS 

Electromagnetic theory shows that the light ex- 
tinction coefficient / u is uniquely determined by the 
thermodynamic properties of the flow. Its space-time 
correlation might therefore be used to study the tur- 
bulent fluctuations of thermodynamic flow properties. 

Consider the experimental arrangement in Figure 
i. Two narrow beams of white light traverse a jet in 



FIGURE 1. OPTICAL MEASUREMENT OF THREE 
DIMENSIONAL WAVE NUMBER COMPONENTS 

the y and z directions. Using an optical filter in front 
of the photo-detectors, we limit the received radiative 
power to the wavelength interval AA centered around 
the wave length A. The light extinction coefficient fi 
then accounts for the power loss through the interval 
Ax centered around the point 5fand is defined by 

. A . Filtered power loss along Ax 

fi(x, t, A, AA; = — r t, r- 6 — ” •( 1) 

p Filtered energy flux . Ax • A A ' 

* Research Physicist, Illinois Institute of Technology 
Research Institute. 



The tui’bulence causes temporal fluctuations 

ix ’ (jT; t, X, AX) - ix- <ix> (2) 

which lead to a fluctuation of the received spectral ra- 
diance (power /solid angle, area and wavelength in- 
terval) 


I (t) = < I > + i (t) . (3) 

y y y 

The intensify fluctuations of the y and z beams are now 
correlated. The resulting normalized temporal cross 
correlation function G is then related to the light ex- 
tinction coefficients ij . ’ by an integral which can be de- 
rived by calculating the time histories of I (t) and 

I (t) through a quadrature of equation ( 1) . ^ 
z 


GU , T) r = 


< i (t) i (t + r)> 

z y 

< I > . < I > 

y z 


+ « 

= / / (x, y, z + £, t) + 4, y + rj, z, t + t) > drj d£ . 

_ 00 


(4) 


This area integration is extremely useful in many ap- 
plications as will now be demonstrated by interpreting 
equation (4) in tex’msof three -dimensional wave num- 
ber components. Consider the fluctuations around a 
point jf, which is defined by the intersection of the up- 
stream beam and the line of minimum beam separa- 
tion. The space -time correlation function around this 
point, 


(5) 


< p' (x, y, r, t) p' (x + {. y + 7), z + t. t + r; > , 


could then be used to describe the correlated fluctua- 
tions in its vicinity through a sy stem of standing waves. 
«Ve are interested in the wave fronts which are parallel 
to the beams. The contributing waves are character- 
ized by the spatial and temporal wave numbers (fre- 
quencies) k, (,t. The mean square amplitude ("en- 
ergy") E of these wave number components is then 
given by [i] 

V k> “V* (6 > 

+cO 

+C0 /* c ». 

1 / f e _i(kx + WT) ( II R iM ,T) x d<1 d * dT * 

8 ii 3 JJ -oo M x 

-00 

The inner integral resembles the measurable quantity 
G, equation (4) , very much. The functional similarity 


is now made complete by assuming that tjiere is. at 
least one line inside the wave front, along which the 
turbulent fluctuations are spatially homogeneous. 
Without loss of generality, this line could be the £ 
axis, since a rotation around the wave normal will not 
change equations (4) and (6). By malting thisassump- 
tion, equation (5) becomes 

<p'<x,y,zH,t)ju'(x + £,y + t;,z,t+T)> 

= </i , (x,y,z,t)p'(x + 4,y + ij, z-f,,t+r)>. (7) 

Substitution of equation (7) into equation (4) and trans- 
formation of £ to -£ shows that the measurable signal 
is equal to the area integral over the space-time cor- 
relation function 

+00 

G(«. t) ? « II R ix ( ?’ T) x dT?d? - (8) 

-oo 


The wave number component E can thus b measured 
directly by repeating the cross correlatio». of light in- 
tensities for various beam separations |. Substituting 
equation (6) gives the final result: 

+oo 

E « (k,W> x = TT // G( £> T ) e " i(k * +a,r)d £ dr * < 9 > 

Obviously, the measurement of three-dimensional 
wave number components is not restricted to the 
3treamwise direction. Exactly the same result holds 
for any other wave front to which the beams can be 
aligned parallel. Furthermore, the cross correlation 
c ould be measured for other light wavelengths X and 
wavelength resolutions X/AX by recording more than 
two phototubes either simultaneously or in repeated 
runs. In this way, one obtains an absorption and/or 
emission spectrum of the wave number component 
E (k, io) , which could then be resolved into the wanted 
thermodynamic property. 

H. MACH WAVE SOUND EMISSION 

To demonstrate the potentialities of the optical 
cross correlation, we now discuss a well known ex- 
ample, which has so far defied a direct me? surement. 
This is Mach wave sound emission. Consider the mean 
square pressure fluctuations <p 2 (xj) > in the far field 
of a supersonic iet as shown in Figure 2. These pres- 
sure fluctuations receive a contribution from Mach 
wave sound sources at position x* on the jet, which has 
been given by Ffows Williams [2], In our notation, 
his superposition of Mach wave sound is given by a 
volume integral 



(x f -x) 2 (z f -z) 2 




2 

F(x) dx* 


( 10 ) 
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FIGURE 2. MACH WAVE SOUND EMISSION 

ov'jr the volume integrated space-time correlation 
function of the density time derivatives. 


Fix) 




>t»r. on 


Here a Q denotes the far field velocity of sound. The 
problem is to measure the combined source strength 
or forcing function F. Obviously, the beams should 
be aligned parallel to the front of the Machwaves; that 
is, we turn the coordinate system ( £, 77 , f) until the £ 
axis coincides with direction of Mach wave emission. 
In this case the vectors "xj -It, andTf -x -"fare al- 
most parallel to the £ axis, as shown in Figure 2, and 
we may simplify the retarded times in equation ( 11 ) 
by the approximation 

|^-x | - | lf f -X -f | = (12) 


For stationary time series , the absolute values of the 
retarded times are unimportant; only their difference 
4 ;/a 0 affects the correlation. Furthermore, the cor- 
relation of density derivatives can be expressed by the 
second derivative of the density correlations. 



9* 

- 972 <P’ d) P' (t + T)> . 


(13) 


The Mach wave sound might therefore be expressed in 
terms of the directly accessible density cross corre- 
lation function 


+» 8 j ♦« ^ 

F(x) - J — JJ <p' (x,t)p'(x + £, t + r (|)>dqd{) <U . - (M) 
-« 8t 


where r (£) « The inner integral is identical 

with the area -integrated correlation function of equa- 
tion ( 8 ) , except that the extinction coefficient fluctua- 
tions p' have been replaced by th density fluctuations 
p'. Inmost flows with frozen composition the specti’al 
region of the detected light may be adjusted to make 
these fluctuations directly proportional : 

p' = C(A, AA) p’ . (15) 

In this case, the inner integral is directly proportional 
to the normalized optical correlation function G. The 
Mach wave sound sources follow from substituting 
equations (15) and ( 8 ) into equation (14): 

+0 ° 8 2 

F(x) = C 2 J — z G(|, r(*))~d* . ( 16 ) 


The argument ({, r ( 4 ) ) describes the temporal cor- 
relation between the space-fixed position If and one 
movingwave front. Foralmost frozen patterns of tur- 
bulence, the temporal changes are concentrated at tl. ' 
space-fixed observer and might be approximated with 
the flashlight image (t= 0 ) of all waves along the same 
normal. All of these waves could be recorded simul- 
taneously by aligning many beams parallel to the tj 
axis, or, by using one thick "y" beam. Therefore, a 
quick look at the exact line integral of equation ( to) 
might be obtained from a "on? shot estimate" emp’oy- 
ing a beam , with a slit shaped cross section. 

.'2 

F(x) at C’ 2 -~j C* (t = 0 ) . (17) 

d'i ‘ 

The asterisk denotes that theslitted section of one 
beam must be enlarged to cover the whole correlation 
volume around the point x*. 
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VII. STRUCTURAL DYNAMICS 




ELASTIC STABILITY OF A SLENDER BAR UNDER DYNAMIC LOADS 

By 

Frank C . Liu 



SUMMARY 


This research is concerned with the elastic sta- 
bility problems of a space rocket under nonconserva- 
tive, dynamic loads. The mathematical 'analysis is 
made on a simple model, a uniform bar with free-free 
ends with one end subjected to a time-varying thrust 
load. This work consists of three parts: The first 
presents the critical thrust of the bar together with an 
outline of various other buckling criteria under differ- 
ent assumptions; the second shows the elastic stability 
of a spacecraft due to the coupling of longitudinal and 
transverse vibrations; and the third part deals with the 
critical thrust of the bar subjected to highspeed aero- 
dynamic load. The variation of frequency of transverse 
vibration due to the thrust and flight velocity is also 
illustrated. \ 


solid rod is equal to the speed of sound in the rod; in 
general, the frequency of the stress wave traveling 
back and forth in a solid metal rod is much greater 
than the natural frequency of transverse vibration of 
the rod, and hence, the coupling effect is negligible. 
However, a rocket usually has great mass density with 
a small effective load-carrying cross-sectional area. 
Its longitudnial frequency may be of the same order of 
magnitude as its transverse vibration frequency. 
Thus, the elastic coupling becomes significant. 

Second , we are concerned with the elastic stability 
of the space vehicle under the combined action of a 
thrust and an aerodynamic load. For a space rocket 
of large cross-sectional area, the aerodynamic force 
and damping effect caused by its transverse oscillation 
in supersonic air flow may have some significance to 
its dynamic stability. Such information may be help- 
ful to structural design engineers. 


I. INTRODUCTION 


II. THE CRITICAL THRUST 


The means of finding the static buckling load of a 
slender column is the well known Euler's method. If 
the external forces have a potential, they are called 
conservative forces, in which case this method is ap- 
plicable [i]. The theory of elastic stability of non- 
conservative problems is based on the investigation of 
the small oscillation of the system about its equilibrium 
position. Thorough discussion and many practical ex- 
amples of the so-called dynamic method may be found 
in V. V, Bolotin's books listed as References 1 and 2, 
This approach is now used to deal with elastic stability 
problems of space vehicles. 

The s. ace vehicle is treated as a uniform, slender 
cylindrical bar, and such that the shell effect, local 
budding, and panel flutter are disregarded. We are 
mainly interested in the elastic stability of the vehicle 
under dynamic loads. This bar has both ends free and 
is subjected to a thrust force at the tail end along the 
direction of the body axis; i.e, , the rocket engine 
gimbal is clamped. Thus, this load is nonconserva- 
tive, because it vibrates together with the bar. Two 
dynamic stability problems are presented here. 

First, we consider the dynamic coupling of longi- 
tudinal vibration and transverse vibration of the ve- 
hicle, The velocity of stress waves traveling in a 


Presented is the analysis of the dynamic approach 
on elastic stability of a slender bar while the well- 
known static approach is also given for comparison. 

A. THE DYNAMIC ELASTIC STABILITY 


When the load is not acting along a fixed di- 
rection due to the motion of the bar, equilibrium of the 
bar cannot be obtained without taking into account the 
inertia force. In other words, the dynamic elastic 
stability of the bar is governed by the differential 
equation of transverse vibration f 2 ] , 


aV 


aw , 


EI 8? + (N + m 


a 2 w 
at 2 


= o 


(i) 


where N = P 0 (1 - x/L) is the axial force along the 
bar, EI is the bending stiffness, m is the mass per 
unit length, and Lis the length. Two types of boundary 
conditions are treated. 

Case 1. The motion of two ends constrained 
along a fixed direction: As shown in Figui’e ie, the 

solution of equation i gives the critical load^ 3 ] 

P* = 2P 
o 
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where P = 7r 2 EI/L 2 is the buckling of a bar with hinged 
ends (Fig. ia). 




FIGURE i. THE CRITICAL THRUST 

Case 2 , Bar with free -free ends: If the bar 
is treated as free -free, we may use the bending mode 
of a free-free bar instead of a sinusoidal mode shape. 
The solution of e nation 1 yields t 3 ] 


P 0 = 8.23P. 


B. THE STATIC ELASTIC STABILITY 

If the lateral vibration of the bar is disre- 
garded, the deflection of the bar is due to the inertia 
load of uniform acceleration. Consequently, the de- 
flection does not depend on time, and the differential 
equation to be dealt with is t 4 ] 

where q = m(P 0 /mL) is the inertia lead. The critical 
thrusts of the bar for various boundary conditions are 
given as follows: 

Case 1. w(0) = w’(0) = 0, as shown in Fig- 
ure lc. The critical thrust given by Reference 4, p. 
116, is 

P* = 0.794P . 
o 

Case 2 . w"(0) =w"(L) =0, as shown in Fig- 
ure Id. Reference 5 gives 

P* = 2.601P . 
o 

Case 3. w(0) =w(L) = 0, as shown in Fig- 
ure lb. From Reference 4, p. 122, which uses en- 
ergy method , we find 


P = 1.88P 
o 


IH. ELASTIC STABILITY DUE TO COUPLING OF 
LONGITUDINAL AND TRANSVERSE VIBRATIONS 


In Section II N of equation ( 1) is treated as the 
average axial force which does not include the longi- 
tudinal stress wave in the bar. To take into consid- 
eration the longitudinal vibration given by the equation 


EA 3x2 m 8t2 


we have 


8u 


N(x,t) = EA^ , 


with 

N ( 0 , t) =-P G (t) 


and N(L,t) = 0, 


(3) 


(4) 


( 5 ) 


where EA is the extensional stiffness of the bar. The 
above relationships and equation (i) indicate that both 
u(x,t) and w(x,t) depend on the thrust P Q (t) , Gen- 
erally, a rocket engine requires 0, 5 to i. 5 seconds to 
build up to its full thrust power; it then keeps a fairly 
constant magnitude at its full power. Let us repre- 
sent the thrust buildup by a simple two-parameter 
curve in the form 


p 0 (t) = p 0 [(1 + & n/t 0 " 


= p 


( 6 ) 


t a t 


where t Q denotes the thrust buildup time , and the con- 
stant a has the value from 0 to 2, as shown in Figure 

2 . 

p o (t> 



FIGURE 2. THRUST BUILDUP CURVES 

Using Laplace transform and Galerkin's tech- 
nique to solve u(x,t) andN(x,t) from equations (3) and 
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(5), then w(x,t) from equation (1>. we obtain the 
waU^OMwa Mathein equation in the generalized coor- 
dinaXe lilt) 


♦ c?d 


B.23P 


- 2 lit BonQ t t)f t 


where fli and flj are fundamental natural frequencies 
of tranareree and longitudinal vihratlona, respectively, 
and fn la a parametric constant. From equation (7 ) 
the stability boundary of flj/O t vs. the thrust factor 
can be plotted. A typical Illustration is shown In Fig- 
ure 3 for a = 1 and r 0 = «jt = 1, 3, and 6. 


(1) the spacecraft Is treated a» a'hollow, slender' 
cylindrical bar*, 

(3) the shell effect and local buckling are diere- 
garded; i.e. , the. cross section is unde formed; and 

(3) the linearised piston theory is applicable. 

Thus, we may put the aerodynamic pressure and q(x,t) 
respectively, 

* p - , ®* r 

P * ” — ( - V — ) q(x,t) =-2 j pco a9 Bd0 






FIGURE 3. THE ELASTIC STABILITY BOUNDARIES 


Figure 3 shows thatthe Instability region Is larger 
for smaller thrust buildup time t 0 and that the width of 
the Instability band increases with the increase of tike 
thrust P . 


DYNAMIC ELASTIC STABIUTY UNDER 
AERODYNAMIC LOAD 


Now, the term of aerodynamic loadq(x,t) exerted 
ion the spacecraft during high speed flight is added on 
the right-hand aide of equation (i). A simplified ap- 
proach has been used to derive this function under the 
following assumptions: 


and let 

eo 

% ■ F v ( (S£)f 1 (t)cosfl 
i=l 1 1 

where w r is the redial deformation of the shell struc- 
ture, U Is the velocity of flight, y is the air constant, 
P« *« fba atmospheric pressure, c ie the velooify of 
sound in air, vj(x) is the ith freebee bending mode, 
and fj < t) is the 1th generalized coordinate. Solving the 
differential equation by dalerkin’a method and then 
applying Routh-Hurwitz criteria to the system of equa- 
tions of the generalized coor dinate s,' we can c a lc ula te 
the critical thrust as function of flight velocity. An 
illustration of the effect of aerodynamic load on criti- 
cal thrust is shown in Figure 4. 




‘ .V. . i !* i-U‘ ‘ 
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FIGURE 4. CRITICAL THRUST VS. FLIGHT SPEED 

The model used for numerical illustrations in 
Figure 4, as well as Figures 5 and 6, has the follow- 
ing structural data: 

R = 15 ft (4. 57 m) 

L = 350 ft (107 m) 

El = 50x1 0 12 in. 2 -lb ( 1. 46xi0 10 m 2 -Kg) 

mgL =? 6x10® lb (2. 72x10® Kg) 

■ All these curves are calculated by using the first two 
generalized coordinates, fj and ig. 

From Figure 4 we observe that P* = il*23P; this 
; is considerably higher than the value°8. 23P given by 
! given by using fj alone. The critical thrust decreases 
. with the increasing of flight velocity, and is not re- 
: duced significantly by the structural damping. 


A. VARIATION OF FREQUENCY WITH THE 
THRUST 


The frequency of transverse vibration of the 
rocket under the action of a constant thrust varies 
with the magnitude of the thrust as shown in equation 
. (7), The frequency ratio, frequency with thrust P 0 
to frequency with zero thrust, vs. the thrust P 0 is 
plotted for various aerocfynamic parameter 0 as shown 
in Figure 5. The lower curve, which does not de- 
pend on 0 , is obtained when f t is used alone. The 
value of 0 plotted against flight velocity is shown in 
the upper right-hand comer. 



FIGURE 5. VARIATION OF FREQUENCY VS. 
THRUST 

In Figure 6 the variation of frequency is plotted 
vs. 0 with the thrust factor K(= P/P) as a parame- 
ter. 



FIGURE 6. VARIATION OF FREQUENCY VS. 0 
V. CONCLUSIONS 

To sum it up, the following findings may be sig- 
nificant: 

1. The critical thrust of a free-free bar under an 
axial thrust is four timer greater than that of a sup- 
ported-supported bar under a conservative force. 

2. The region of elastic stability of a rocket de- 
creases with the thrust buildup time. 

3. The critical thrust decreases rapidly at high 
supersonic speed and low flight altitude. 

4. Structural damping does not have significant 
effect on the critical thrust. 
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OPTIMIZATION OF SLOSH BAFFLE GEOMETRY 
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SUMMARY 



This paper describes a method for determining 
the optimum baffle configuration for suppressing ex- 
cessive propellant sloshing in the tanks of a large launch 
vehicle. The optimum geometry for a baffle system 
is considered to be the one with the least mass which 
is capable of providing the necessary damping. The 
analysis considers acylindrical tank and a baffle sys- 
tem made up of several flat ring baffles moimted hori- 
zontally* The basic assumptions are discussed and 
derivations givenof the equations necessary to predict 
the optimum baffle configuration consistent with ade- 
quate strength. Numerica 1 results are presented 
graphically, and some general conclusions are drawn. 


D Vertical spacing between baffles 

Fl,F 2 ,F 3 Functions defined by equations (ll,i2,5) 

h Length of portion of tank where damping is 

required 

H Total length of tank 

K Parameter which describes the critical 

baffle depth 

1 Distance from center of tank to center of 

pressure measured radially 

P Load applied to baffle due to pressure 


I. INTRODUCTION 
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The phenomenon of propellant sloshing is particu- 
larly important to the design of large liquid-fueled 
launch vehicles due to the interaction of fuel sloshing 
with the vehicle’s control system. The simplest and 
most widely used method for reducing the effects of 
propellant sloshing is by means of a system of baffles. 
Several baffle configurations have been employed; one 
of the. most common has been the flat ring baffle ex- 
tending inward from the tank wall. Since high vehicle 
performance requires that the mass of the vehicle 
structure be kept at a minimum through efficient de- 
sign, it is highly desirable to use the most efficient 
baffle design for the particular propellant tank. 

This analysis outlines the basic characteristics of 
an efficient baffle system and gives insight into the 
manner in which the weight varies with change in con- 
figuration of the syst^' 


I LIST OF SYMBOLS 

i 

j 

Symbol Definition 

, A Area of single baffle 


n 


N 

P 


r 


t 

v 

w 

wt. 

a 


Mass of fluid striking the baffle 

Mass rate of flow of fluid striking the 
baffle 

Number of submerged baffles contributing 
to damping 

Total number of baffles in the system 
Pressure on baifle 

Radius from tank ..enterline tu a point on 
the baffle 

Design stress in tension of baffle material 
Baffle thickness 

Velocity of fluid at impact with baffle 
Baifle width 
Baffle weight 

Longitudinal acceleration of tank 

Total fluid damping provided by n sub- 
merged rings 




Tank radius 

Depth of nth baffle below equilibrium 
surface 


t Maximum amplitude of fluid oscillation at 

v tank Wiill 

cj Natural frequency of fluid oscillation 


\ 

t 

i 

j 

\ 


9 r > 


I 



Mass density of the fluid 


B. DERIVATION OF OPTIMIZING EQUATIONS 


P L 

pg Mass density of the baffle material 

v Poisson's ratio 

II. ANALYSIS 

A. BASIC ASSUMPTIONS 

The basic tank configuration is shown in Fig- 
ure 1 below. The tank is cylindrical with a flat bottom 
and the baffle system is made up of a series of evenly 
spaced, solid rings mounted horizontal 'y. These as- 
sumptions are based on the conclusion made by Langner 
[1] that for moderate damping a solid ring baffle sys- 
tem has minimum weight. Li arriving at these con- 
clusions, Langner conducted a preliminary analysis of 
ring and partition baffles, both perforated and solid. 
He also concluded that, with the exception of very thin 
baffles where the thickness is limited by manufacturing 
tolerances, there are no distinct savings in weight to 
be gained by perforating ring baffles. The analysis 
provides for an arbitrary region at the bottom of the 
tank where damping is not required because of the 
small amount of liquid remaining in the tank. 



FIGURE i. TANK AND BAFFLE SYSTEM 
GEOMETRY 


1. Damping Law . A damping law for a flat 
ring baffle mounted in a cylindrical tank was first pro- 
posed by Miles [2] ; comparison with extensive meas- 
urements [3, 4, 5) has confirmed that this relation- 
ship gives a resonable approximation of the actual 
damping. Bauer [6] expanded Miles' equation to ac- 
count for more than one baffle. In this analysis the 
contribution to the total damping by exposed or par- 
tially exposed baffles is taken into account through the 
use of an experimentally determined parameter, K, 
including the combined effects of exposed atid unex- 
posed baffles. K is defined below in terms of baffle 
depth and spacing. The depth, d ls of the first baffle 
beneath the undisturbed free surface when the tank is 
filled is indicated in Figure 2; dj corresponds to the 
fluid levol at which the damping provided by the first 
ring is equal to the minimum provided by the first and 



DAMPING RATIO, y ( 


FIGURE 2. TYPICAL VARIATION OF DAMPING 
RATIO WITH FLUID LEVEL 

second rings combined. In this manner, it is possible 
to insure that the damping provided never falls below 
a prescribed value. A relationship between dj and D 
can be written < 

d,/D= K, 

where K is estimated from experimental data. Then * 
the depth, d^of the nth ring below the surface is given 
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by 


d n = KD + (n-i) D. 


the following procedure can be used to estimate this 
pressure distribution. If the fluid intersecting the 
( i) baffle is completely stopped, 


P = m v. 

The total damping provided by the system of n rings is 

given by Bauer's expansion of Miles' equation. where 




n/ i^/a e -4 * 6 



( 2 ) 


By substitution for d n , equation (2) can be rewritten 
as 


Y T 



2.83*7 £, /a 
w 


= f) e-4.6 [KD/a+ (n-l)D/a] 
n=l 


V=fw m = p L AV Then p =~= p L (t,w ) 2 . 

From f4 ] 


co = J ~ e tanh (e h/a) 

where for the first mode e = 1. 84119 and for the case 
of h/a ^ 1.0, tanh (e h/a) ^ 1.0. Now by assuming 
that the liquid surface remains planar, the displace- 
ment from the equilibrium position, £ , at any point is 
related to the displacement, £ w> at the wall by 

£ = r / a S w • 

The pressure distribution on the baffle is given by 


The next step is to determine the dimensionless spac- 
ing, D/a, as a function of the dimensionless width, 
w/a, for given values of the constants K, n and 

^/a ; however, because of the form of equation 

(3), this is difficult to do analytically. By resorting 
to a curve-fit procedure, a fifth order polynomial for 

D/a in terms of w/a and y rj y*7 ^/a. can be found: 


D/a = a 0 + aj F 3 + a 2 F§ + a 3 Fjj + a 4 f| + a 5 F|, (4) 


where 



2. 83 *7 £, /a 
w 



-f> 


1-3/2 


(5) 


The coefficients for the terms in equation (4) are de- 
termined from the method of least squares and are 
functions of the values selected for the constants K t 
and n. 


2. Pressure Distribution on Baffles . Data 
from large scale slosh tests [7] indicate that the 
maximum pressure differentials on ring baffles are 
experienced when the equilibrium fluid level coincides 
with the baffle position. It is thought that this pressure 
is the resultof the oscillating fluid striking the baffles 
from beneath with its maximum kinetic energy. The 
maxi mum pressures for this case are found to lie 
along the tank diameter parallel to the direction of 
excitation. Assuming that only the first mode of fluid 
oscillation contributes to the pressure on the baffles, 


P T <* I £ w \ 

P= 1.7506 |-~J (r) 2 . (6) 

Figure 3 shows some of the data from Reference 7 
compared to the theoretical values obtained from 



FIGURE 3. BAFFLE PRESSURE VS. AMPLITUDE 
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equation (6) . The predicted pressure is of the cor- 
rect magnitude, but the variation of pressure with fluid 
amplitude has not been accurately described. This can 
probably be traced to the assumption of a planar sur- 
face, since after the wave has struck the baffle, its 
surface would obviously be broken up and made more 
complex. The degree to which the experimental data 
can be grouped by dividing by the radius squared sug-, 
gests that the shape of the distribution on the baffle is 
essentially correct. 

A the or, for predicting the pressure distribution 
as a result of continuous flow past a ring baffle has 
been presented by Liu [9]. While this method seems 
to work well for the case of the baffle being completely 
submerged, Junkins [7] showed that for the case be- 
ing discussed here the pressure prediction fell well 
below measured data. Therefore, for the purpose of 
gaining insight into the optimization problem, equa- 
tion (6) will be considered adequate. 

3. Loading due to Pressure. The pressure 
distribution given by equation (6) applies along the di- 
ameter parallel to the direction of excitation as indi- 
by Figure 4. The relative magnitudes of the pressures 



FIGURE 4. PRESSURE DISTRIBUTION ON A RING 
BAFFLE 

in the distribution are greater than those on any other 
diameter since the fluid amplitude at the wall is a 
maximum. To apply the existing equations for stresses 
in flat rings [8], it was necessary to assume that this 
distribution was made uniformly on all diameters. 
This assumption is somewhat conservative but is a 
reasonable approximation of narrow rings. The total 
load due to this pressure distribution can be found by 
integrating over the entire baffle. 


a 

P = 2 7 r J p r dr = 



where 


rj = a-w. 


Then 

p = 2.748 a 3 p T a 

Xj 


4 

[1 - (1 -w/a)] . (8) 


The center of pressure can be found in like manner 
by integrating the moment over a unit arc length. 



0.800 


fl- ( 1 -w /a ) 1 . 
fi- (i-w/a) 4 ] 


(9) 


By assuming that the baffle derives its principal 
strength from resistance to bending, it is possible to 
make use of the equations for bending stress as ap- 
plied to the equivalent system shown in Figure 5 . After 


i 



i 


FIGURE 5. EQUIVALENT MANNER OF LOADING 
(OUTER EDGE FIXED) 

substitution of equations (8) and (9) , this appears in 
the form 

a %T a <V a > 2 

s, - — ~ — 0.65616 Fjd -F 2 ) J, (10) 


i 
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where 


(- + 1 ) 

v 

1 

V 

i i 

_C£ 

a 


(log 


1 /a 
cp 


i* 


F 2 = 


- 1 ] [i- (1 -w/a ) 4 ] 

(£- 1 ) -( 7 +I) ( 1 -w/a ) 2 
( — -!) + (— + 1 ) ( 1 -w/a ) 2 


_ V 


examined using the methods developed above are dis- 
cussed briefly. Since there is some question as to the 
validity of the damping law when applied to very wide 
baffles, only values of w/a up to 0. 20 have been in- 
cluded. However, this includes the range of practical 
interest since nothing seems to indicate any advantages 
for very wide ring baffles. The constant K was as- 
signed a value of 0. 55 after an examination of the ex- 
perimental data of Reference 5 . Poisson's ratio 
was assumed to be equal to i/3 since this corresponds 
roughly to most metals that might be considered. 


If equation ( 10 ) is solved for t, the result is given by 



a 


0 . 

s t L 


65616 F^l + F 2 ) .(13) 


']* 


✓The total weight of the baffle system is 


Figure 6 shows the typical variation of the weight 
parameter with w/a. The shaded area indicates the 
range of values for w/a which cannot produce the re- 
quired damning regardless of the number of such baffles 


Wt. = N p a t [t a 2 - it (a -w) 2 ], (14) 

X) 

However, the number of baffles, N, can be written in 
terms of their spacing as 


N = h/ti . 


By substituting this and the expression for t given by 
equation ( 13) , it is possible to write a dimensionless 
weight parameter of the form 
W t. 

7 — me 

a p B a s (h/a) ( t w /a ) v ' P L s t 


» 0.810 TT 


1 

2 


[Fid + F*)] 


[1 - (1 -w/a)‘] 


D/a 


(15) 


where D/a is given by equation (4) . Since this weight 
parameter contains all the arbitrary quantities perti- 
nent to the baffle material and the fluid, the configu- 
ration for which the parameters are a minimum will be 
the optimum baffle configuration regardless of the ma- 
terial or fluid to be used. The only parameter ap- 
pearing on the right side of equation (15) which is not 
a function of the geometry of the system is 


Y T 


/a 


in. DISCUSSION OF RESULTS 


The results of a number of cases which have been 



FIGURE 6 . NONDIMENSIONAL WEIGHT PARAM- 
ETER VS. NONDIMENSIONAL WIDTH 


used. A value of 0 . 10 was selected for the parameter 
y T 

which contains both the desired damping and 

the maximum fluid amplitude. In addition , only the first 
three submerged baffles were assumed to contribute to 
the damping. The optimum baffle geometry produces 
the minimum value of the weight parameter as indi- 
cated in Figure 6 . By constructing several such 
curves and by selecting from each the minimum value 
of the weight parameter, the optimum geometry may 

Y T 

be presented as a function of the parameter ~ . 

^£ w /a 
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The variation of the optimum value of nondimensional 
width and nondimensional spacing with - — Is 

7 a ' 

shown in Figures 7 and 8, respectively. Figure 7 in- 
dicates that as the parameter increases so 

■TQP 



FIGURE 7. VARIATION OF OPTIMUM NONDI- 

y 

MENSIONAL WIDTH WITH T - 

does the optimum value ofw/a. Figure 8, however, 

seems to indicate that the optimum value of D/a is a 

y 
r t 

constant for all the values of considered. 
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FIGURE 8. VARIATION OF OPTIMUM NONDI- 

y 

MENSIONAL SPACING WITH T 
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SATELLITE SHADOW TIME 


-223 *4 ! 

P, E. Dreher 


SUMMARY 

/ ' < : 

Analytical methods have been developed to deter- 
mine the time spent in shadow by an earth -orbiting 
satellite. The length of the shadow time is a function 
of the particular orbit under consideration (i.e. , semi- 
major axis and eccentricity) and of the orientation of 
the orbit with respect to the sun. The technique as de- 
veloped was programmed on the GE-225 and a general 
parameter study for two different orbits was per- 
formed. The umbra transit times are as great as 12 
hours for certain orbit orientations and less than 0,5 
hour for other orientations. The penumbra transit 
time is generally less than 30 percent of the umbra 
transit time. 

The effects of the refraction of the sun’s rays by 
the earth’s atmosphere on these transit times can be 
significant. Since the refraction method has not been 
incorporated in the computer program, only one case 
was evaluated. For this particular case the umbra 
transit time was reduced approximately 300 percent. 

t 

V 

I. INTRODUCTION 


Presented are some of the results of a study per- 
formed by Northrop Space Laboratories, Huntsville 
Department, for the Operations Studies Branch, Aero- 
Astrodynamics Laboratory, under Contracts NAS8- 
11111 and NAS8-20082, ’’Satellite Shadow Time for 
Earth Orbital Flights,” [1], presents some of these 
results. 

For an earth satellite dependent upon onboard so- 
lar cells for ene jy, the length of time spent in the 
earth’s shadow duungeach orbital revolution is criti- 
cal. Consequently, prior knowledge of the shadow 
history for the orbit under consideration is necessary 
in the design of the onboard energy supply. 

An analytical method has been developed for de- 
termining the time an earth-orbiting satellite will re- 
main in the earth's umbra and penumbra while traveling 
in a specified elliptical orbit. The analytical equation 
for the orbit ellipse under consideration is solved si- 
multaneously with the equation for the shadow cone to 
obtain the intersection points. The eccentl'ic anomaly 


corresponding to each shadow intersection point is used 
in Kepler's equation to determine the shadow transit 
time for the given ellipse. 

The effect of atmospheric refraction on this sha- 
dow time has been investigated using Fermat's optical 
path length principle. The umbra transit time is de- 
termined considering the satellite to traverse the um- 
bra formed by the refracted rays rather than the 
straight rays. 

II. ANALYTIC FORMULATION OF SHADOW 
TRANSIT TIME 


For this analysis the following assumptions were 
made: 

l c The earth is a perfect sphere of radius r^, 

2. All perturbational forces are neglected. 

3. The sun 12 a perfect sphere of radius r^. 

4, The earth-sun separation is a constant 
distance dg, 

5, The relative motion between the sun and 
earth during shadow transit time is 
negligible. 

The earth’s umbra in space can be represented 
by a circular cone as shown in Figure 1. The trace 
of the circular cone on the orbital plane will be a conic 
section, the exact form depending on the relative 
orientation of the orbital and ecliptic planes. To have 
transit time in umbra, the trace of the circular cone 
on the orbital plane must intersect the orbit ellipse 
in two points which are located on the umbra boundary. 

To determine the intersection points of the circu- 
lar cone trace and the orbit ellipse, the analytical 
equations for each conic section are solved simultane- 
ously, Such a simultaneous solution yields all inter- 
section points of the two complete conic sections. The 
two intersection points which lie on the umbra bound- 
ary define the shadow time history . 


It is convenient to define two coordinate systems, 
one in which the circular cone can be easily described 
by an analytical equation and another in which the el- 
liptical orbit can be similarly described. Two such 

systems are the {X Y z}and the {X Y Z } systems 

L L L 
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FIGURE 1. CONE REPRESENTATION OF UMBRA 

shown in Figures 2 and 3. The {X Y z} system is a 
rotating geocentric ecliptic coordinate system with the 
X-axis pointing directly away from the sun, the Y-axis 
pointing to the ecliptic north, and the XZ -plane being 


y 



y, the angle between the ascending orbital node 
and Rp, the vector through perigee, measured in the 
orbital plane. 

The standard equation for the circular cone in the 
{X Y Z} system is 


( X-a t ) i, 2 _ _Y^ 

aj 2 b i 2 



where 


ai 


r E d 

r s " r : 


and 


The standard equation for an ellipse in the {X E Y E ?.£) 
System is 


(X £ + a2e) 2 


2 

»2 



b 2 


where 


- 1 , 


FIGURE 2. {X Y Z} COORDINATE SYSTEM 



FIGURE 3. {X E Y fi Zg} COORDINATE SYSTEM 

the ecliptic plane. The {Xe Ye Z e } system is a ro- 
tating geocentric orbital coordinate system with the 
Zg-axis pointing to perigee, the Zg-axis pointing a - 
long the semi-latus rectum to a true anomaly of 270° , 
the X E Z E -plane being the orbital plane, and the Yg- 
axis completing the right -hand system. The orthogonal 
transformation from the (X Y Z} system tothe fX E Y E 
Ze } system is developed by means of three successive 
rotations through the angles: 

p , the angle between the sun and the ascending 
orbital node (the intersection of the ecliptic and the 
orbitr’ t >; m 3 ) measured in the ecliptic plane 

i, the ino ination between the ecliptic and the or- 

bital plane 


a 2 


b 2 


e 


R 


i«vy 

^r x a 

a p 


apogee distance 
perigee distance. 


The cone equation can be transformed to the coordinate 
variables of the {Xg Y E Z E } system. The resulting 
transformed equation in the orbital plane, i.e. , the 
X E Z E - plane, and the standard equation for an el- 
lipse in the {X E Yg Zg} system are solved simultane- 
ously to yield all intersection points of the cone traoe 
and orbit ellipse. Since the sun is always looated on 
thenegative X-axis (by coordinate sy stem definition) , 
ashadowpoint must have a positive X-oomponent. To 
select the shadow points, it is sufficient to transform 
all the intersection points to the {X Y Z} system and 
eliminate thosewithanegative X-oomponent, The two 
shadow points are designated as (Xgj, 0, Z E1 ) and 
(Xgg, 0, Z E2 ) . The satellite will spend time in um- 
bra for the orbit under consideration if and only if two 
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real shadow points exist. The eccentric anomaly cor- 
responding to each shadow point is 

„ * h- r '\ 

E > ■ cos 


E 2 = cos 
where 


/a 2 - r 2 
[ a 2 e 


* 


- *ei 2 + z ei ! 

r * ■ ^7 + Z E2 2 • 

The transit time from perigee to a shadow point 
is determined from Kepler T s equation 


bf = 


r„ d„ 
E S 




This technique as presented was developed for 
general studies of shadow time histories as a function 
of parameters defining the orbital plane orientation with 
respect to the sun. To study specific trajectories, 
transformations must be established between the 
{X Y Z} system and the systems in which trajectory 
data are available, e.g. , the {X„ Y E Z E } and (X' Y' 
Z"} systems. The {Z ' Y' Z ' j system is a space - 
fixed geocentric ecliptic coordinate system with the 
X'-axis pointing to the vernal equinox, the Z " -axis 
pointing north and the X' Y" -plane being the ecliptic 
plane. The relative orientation of the {X' Y' Z'}, 
(X Y Z} and {X E Y E Z E } systems is shown in Figure 
5. 


T. 

x 


E. - e sin E. 
1 i 



i = i, 2 


where K = earth’s gravitational constant. From these 
times the umbra transit time DTp can be determined: 


DT fi = | T t - T 2 | . 


A similar procedure is followed to determine pe- 
numbra transit time because the earth's penumbra is 
also represented by a circular cone (Fig. 4) : 



FIGURE 4. CONE REPRESENTATION OF 
PENUMBRA 


(X-ap)‘ 


a i 



where 
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r_ d 
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- 2 
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= 0, 



FIGURE 5. RELATIVE ORIENTATION OF {X Y Z} , 
{X'Y'Z'}, AND {X_, Y Z_} COORDINATE 

ili £4 £4 

SYSTEMS 

The inclination, i, can be determined by a spher- 
ical trigonometric consideration of the intersection of 
the equatorial, ecliptic, and orbital planes as shown 
in Figure 6. The inclination of the orbital plane, I, 
the inclination of the equatorial plane to the ecliptic, 
p, and the angle between the line of nodes and the as- 
cending node, L, are known quantities. 

Once the available trajectory data are used to de- 
termine the parameters 0, y, and i, these param- 
eters , along with the specified apogee and perigee 
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distance R and Rp, can be used to compute the sha- 
dow transit time. 



FIGURE 6. EQUATORIAL, ECLIPTIC AND 
ORBITAL PLANES 


Because the inclination, i, generally varies only 
slightly with time, a general parameter study was 
done varying 0 and y and holding i constant. Two dif- 
ferent orbits were considered with the results shown 
in Figures 7 and 8. 



FIGURE 7. GENERAL PARAMETER EFFECTS 
FOR 1/8 LUNAR PERIOD ORBIT 

Figure 7 shows the umbra transit time over the 
range of 0® to 360® for /? and over the range 202® to 
290* for y. The orbit considered was a preliminary 
1/8 lunar period orbit with apogee « 186,661 km and 
perigee = 7073 km. On the launch date, y - 202® , 
and one year later, y * 290®, As seen in Figure 7, 
the umbra transit time varies from 0. 5 hours to 5.6 
hours depending on the value of /? . 



o 90 m m mo 

FIGURE 8. GENERAL PARAMETER EFFECTS 
FOR 1/3 LUNAR PERIOD ORBIT 

Figure 8 shows the umbra transit time for a pre- 
liminary 1/3 lunar period orbit with apogee = 350, 631 
km and perigee = 7410.9 km. A general parameter 
study was done for this orbit over the range of 0° to 
360° for 0 and 90° to 180° for y. The maximum sha- 
dow time occurred for y and /3 = 180; this is logical 
since this placed the sun on f he orbital node and placed 
apogee in shadow (dire* A y behind the earth from the 
sun) . This case represents the maximum shadow time 
possible for a given orbit. This was not the case for 
the 1/8 orbit considered because actual y conditions 
were used, and the minimum value of y wa3 202°. 
This meant that apogee was not directly behind the 
earth from the sun. The maximum shadow time is 
also a function of apogee altitude. A difference in 
apogee altitude from 186,661 to 350,631 km changed 
the maximum shadow time from 5.6 hr to 12.2 hr. 

Figure 9 shows the transit time through umbra and 
the transit time through umbra plus penumbra. The 
time spent in partial shadow (penumbra) is generally 
less than 30 percent of the time spent in umbra. 


TRANSIT TIKE (HR) 



FIGURE 9. COMPARISON OF UMBRA AND 
PENUMBRA TRANSIT TIMES 
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Figure 10 shows the umbra transit time for a 
particular trajectory, the 1/8 lunar period orbit. For 
this case actual trajectory points were used, and the 
values of i, /?, and y, corresponding to these points 
were computed. The maximum shadow time was 1. 82 
hours with the minimum 0. 36 hours. 


SHADOW TIME (HRj 



FIGURE 10. UMBRA TRANSIT TIME FOR 1/8 
LUNAR PERIOD ORBIT 


HI. REFRACTION EFFECT ON TRANSIT TIME 


In the previous discussion, the light rays from the 
sun are assumed to travel in straight -line motion. But 
the rays are refracted as they pass through the earth's 
atmosphere because they travel betweenmedia of dif- 
ferent densities. The resulting path is one which sat- 
isfied Fermat's principle for optical path length. The 
umbra is no longer a circular cone whose XY trace is 
a straight line , but a surface whose XY trace is the 
curve describing the path of a refracted ray, as shown 
in Figure 11. For simplification, this portion of the 
study was confined to the two-dimensional space, i.e. , 
the XY -plane. 


Y 



FIGURE li. REFRACTED RAY FORMATION OF 
UMBRA 


The following assumptions were made: 

1. The earth's atmosphere is formed by an in- 
finite number of successive strata, each of uniform 
density. By passing to the limit in which the number 
of strata becomes infinite and the differences in den- 
sity in two consecutive strata infinitesimal, one has 
the case of the continuously varying density of the at- 
mosphere. 

2. The successive strata are all concentric with 
the center of a spherical earth. 

3. The values of the atmospheric density at the 
surface of the earth [2] are 

p Q = 0,0027 slug/ft 3 . 

4. The value of the refractive index at the surface 
of the earth [2] is 

n. = 1,00029 . 
o 

5. The atmospheric density model can be ex- 
pressed as 

p (r) ~ — . 

e 

6. The relation between the refractive index and 
the atmospheric density is given by Sellmeier's equa- 
tion f 6] 

n 2 (r) = 1 + C p (r) , 


where 



7. The upper limit of the atmosphere is 
h =6419.719 km. 

To determine the refracted ray equation, Fermat's 
principle is used. The path of a light ray is such that 
the integral 

1 = nds 
A 

is an extremum [4] , where 
1 = optical path length 
n = refractive index 
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ds = differentia?, of arc length along the path A B. 
Calculus of variations techniques were used to deter- 
mine the optical path, this path being that of the re- 
fracted ray. 

As seen in Figure 12, the refraction angle | can 
be defined as 

f = 6 -90 # - 0 , 

s o o 


the orbit ellipse equation will give the intersection 
points necessary for the computation of umbra transit 
time. 

Numerical result? were obtained for one trajec- 
tory with R & = 186, 360. 0 km and Rp = 6545. 0 km. In 
this two-dimensional case, apogee was considered to 
be directly opposite the sun (/? = 180°, ? = 180°) so 
that the shadow time neglecting refraction was the 
maximum possible for this orbit. The results were 


where 

<t> = cos 

o 


n 

o 


E 


h * 


Umbra time, neglecting refraction 8. 03 hr 


) 


Umbra time, considering refraction 3. 22 hr (0=1,2') 


*< 



FIGURE 12. MODIFIED PATH OF THE 
REFRACTED RAY 

The angle 6 can be determined from a definite in- 
tegral that is a function of the refraction index, n, the 
atmospheric density model, and the constants n , r^, 
and h. This integral can be solved numerically 8y us- 
ing Simpson's Method or the Modified Gaussian Method, 
The value of the refraction angle, £, was found to be 
38.8 . 

The other variable used in the analysis is 0 , the 
angle made by the limiting ray with the horizontal be- 
fore refraction, hi the case where refraction is ne- 
glected, the limiting ray makes the maximum angle 
with the horizontal as shown in Figure 12, The max- 
imum value of 0 was determined geometrically, and 
found to be 0^^ =15', 

Figure 12 shows that the ray which was refracted 
by the atmosphere proceeds in straight-line motion 
after leaving the earth's atmosphere ( at point C ) . An 
analytical equation describing the resulting straight 
ray in the {X Y} system solved simultaneously with 


Umbra time, considering refraction 2. 67 hr (0=6. 2') 
Umbratime, considering refraction ?. 3 hr (0=11.2') 

Since this part of the analysis has not been incor- 
porated in the computer program, no other cases have 
been studied. It is not known whether consideration 
of atmospheric refraction reduces the umbra transit 
time by such a large amount for other orbit orienta- 
tions. 

These numerical results show that atmospheric 
refraction has a major effect on the formation of the 
earth's umbra. This should not be surprising, for 
during a total lunar eclipse the moon remains dimly 
visible because some of the sun's rays are refracted 
by the earth's atmosphere and reach the moon [5], 


IV. REFERENCES 


1. Satellite Shadow Time for Earth Orbital Flights. 
Northrop TM, no. 88, 

2. White, J. F. , ed. : Flight Performance Handbook 
for Powered Flight Operations. John Wiley & 
Sons, Inc. , 1963, 

3. Pars, L. A. : An Introduction to the Calculus of 
Variations. John Wily & Sons, Inc., 1962. 

4. Rossi, R. : Optics. Addison-Wesley Publishing 
Company, Inc, , Reading, 1957, 

5. Struve, O. ; Iynds, B. ; Pillans, H,: Elementary 
Astronomy, Oxford University Press, 1959. 

6. Ditchburn, R. : Light. Interscience Publishers, 
1959. 


107 



i 


BLANK PAGE 




IX. PUBLICATIONS AND PRESENTATIONS 





A. PUBLICATIONS 

TECHNICAL MEMORANDUM X-53187 
January 13, 1965 

CISLUNAR METEOROID IMPACT AND PUNCTURE 
MODELS WITH PREDICTED PEGASUS SATELLITE 
PUNCTURES 

By 

Charles C. Dalton 

George C. Marshall Space Plight Center 
Huntsville, Alabama 


ABSTRACT 


The author 's previous model for the variation of 
meteoroid impact flux with respect to distance from 
the earth in cislunar space is simplified and corrected 
for the increase in the inferred mass of a meteoroid 
due to the new pressure values in the U. S. Standard 
Atmosphere , 1962, A new model for the corresponding 
puncture flux is developed as a function of material 
parameters and thickness by interpreting Explorer 
XVI puncture data in relation to the revised model for 
impact flux and to extrapolated results from the Martin 
Company 's multivariate analysis of laboratory hyper - 
velocity impact data. 

During one year after a Pegasus satellite has been 
deployed in an orbit with perigee and apogee heights of 
505 and 740 kilometers, the number of punctures in the 
16 sheets of the 0,0015-inch sensor, the 34 sheets of 
the 0. 008 -inch sensor, and the 366 sheets of the 0. 016- 
inch sensor are predicted to be approximately 5430, 
41, and 38, respectively, when the extra material in 
the bonded capacitors is considered. The estimated 
masses ofthe meteoroids which can just puncture these 
sensors when they have the mean values of density 
( 0,44 gram per cubic centimeter) , velocity (26.7 kil- 
ometers per second), and impact angle ( tt/ 4 radian) 

are iO -7,84 , i0" 6 * 86 , and 10" 4, 98 gram. 


TECHNICAL MEMORANDUM X-53188 
January 15, 1965 

DYNAMIC TEST RESULTS OF SAD -8 ,9 
By 

Gale R, Ernsberger 


George C. Marshall Space Flight Center 
Huntsville, Alabama 


ABSTRACT 


Dynamic tests were performed on a full scale pro- 
totype of the Saturn SA-8 and 9 flight vehicle to deter- 
mine the body bending and torsional dynamic behavior, 
A full scale prototype was vertically suspended in the 
test tower and excited by shakers mounted at the en- 
gine gimbal planes. The vehicle's dynamic response 
at various applied forcing frequencies was recorded 
by vibration pickups. This report presents the SAD-8 
& 9 dynamic test results in terms of resonant fre- 
quencies, mode shapes, generalized masses, damping 
coefficients and response (amplitude and phase) 
curves. 


TECHNICAL MEMORANDUM X-53193 
January 27, 1965 

(MMC) - S-IV/APOLLO SEPARATION AND VENT- 
ING STUDY FOR SA-9 VEHICLE 

By 

D. O. McNiel and W. H. Harmon 

George C. Marshall Space Flight Center 
Huntsville, Alabama 


ABSTRACT 


This report presents the results of a separation 
(MMC) - S-IV/Apollo and venting study for the pre- 
dicted Saturn SA-9 flight. Results in this report are 
subject to the restrictions that residuals are only 
vented in the gaseous state. There are two primary 
objectives of this study, the first of which is to assure 
adequate separation distance between the two vehicles 
[ (MMS) - S-IVand Apollo boilerplate] and the second 
of which is to determine a confidence level that ex- 
cessive angular rates and angular accelerations will 
not occur during and following (MMC) wing deploy- 
ment, Results indicate that predicted (3 a) angular 
rates and regular accelerations do not exceed the 
limits established by the (MMC) Project Office of 6 
deg/sec resultant angular velocity and , 15 deg/ sec 2 
resultant angular acceleration. The maximum pre- 
dicted (3 d) resultant angular rate occurs after all 
re siduals have been depleted and is 4. 7 deg/ sec , Pre - 
dieted resultant (3 «r) angular acceleration during and 
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following (MMC) wing deployment do not exceed 0,005 
m/sec 2 . The minimum (nominal - RSS) separation 
distance between the (MMC) - S-IVand Apollo boiler- 
plate at initiation of wing deployment is 64, which is a 
sufficient distance to prevent collision of the two ve- 
hicles during this operation. 


TECHNICAL MEMORANDUM X-53203 
January 29, 1965 

SA-9 FINAL PREDICTED TRAJECTORY 
PART I. NOMINAL TRAJECTORY 

By 

Gerald Wittenstein 

George C. Marshall Space Flight Center 
Huntsville, Alabama 


ABSTRACT 


This report presents the final nominal predicted 
trajectory for SA-9. It is the first part of a two part 
presentation of the SA-9 Final Predicted Trajectory, 
Part II will present the dispersion analysis associated 
with the nominal trajectory and will be published as a 
separate report (Part II. Nominal Trajectory and 
Dispersion Analysis) . 

The SA-9 vehicle is to be flown over the Eastern 
Test Range (ETR) witha flight azimuth of 105 degrees 
East-of -North. Also included in Part I is a discussion 
of the vehicle and mission objectives, the trajectory 
shaping and constraints, and a description of the ve- 
hicle configuration. 

Primary objectives for SA-9 are to place the Peg- 
asus (MMC) satellite into an elliptical orbit, to further 
flight test the vehicle system, and to further man -rate 
the Saturn class vehicle. The primary objective of the 
Pegasus satellite is to obtain information on micro- 
meteoroids for near-earth orbits. 

The depleted S-IV stage and payload, which in- 
cludes the Pegasus satellite and an Apollo boilerplate, 
are to be inserted into an elliptical, low-earth orbit, 
with perigee and apogee altitudes of 500 km and 750 
km, respectively. The nominal lifetime for this orbit 
is well in excess of one year. 

Information contained in this report may be con- 
sidered applicable for SA-8 and SA-10 flight profiles 
until superseded. 


TECHNICAL MEMORANDUM X-53209 
February 16, 1965 

THE EFFECTS OF INTRODUCING OPTIMIZATION 
INTO A FAMILY OF SATURN IB FIRST STAGE 
TRAJECTORIES 

sy 

John B. Winch and Gary McDaniel 

George C. Marshall Space Flight Center 
Huntsville, Alabama 

ABSTRACT 


This report presents a study of the effects of in- 
troducing lift into a family of Saturn IB first stage 
gravity tilt programs. In this family of trajectories 
the vehicle is constrained to fly a nonlifting trajectory 
profile after a "kick" has been given shortly after 
launch to initiate tilting. Lift Is thus introduced into 
this family of trajectories using p calculus of varia- 
tions technique at 60, 65, and 70 seconds which cor- 
respond to times prior to and after Mach 1 and the 
max q time points. The cutoff weight into orbit for 
these optimized trajectories indicates that an increase 
of as much as 500-600 pounds over the nominal zero- 
lift trajectory may be obtained while the structural 
loadings are still realistic. 


TECHNICAL MEMORANDUM X-53214 
March 8, 1965 

JET PENETRATION INTO A LIQUID 

sy 

H. G. Struck 

George C. Marshall Space Flight Center 
Huntsville, Alabama 


ABSTRACT 


This investigation is concerned with the phenom- 
enon of a gas jet impinging on, and penetrating into, a 
liquid. The study is restricted to the case of a round 
jet of subsonic and supersonic velocities penetrating 
the liquid at right angles. The phenomenon was ana- 
lyzed from two viewpoints. The first, a stagnation 
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pressure analysis, related the depth of the surface 
depression or cavity to the stagnation pressure based 
on the centerline velocity of the jet in the neighbor- 
hood of the surface. The second, a displaced liquid 
analysis which takes the mixing process of the jet into 
account, related the weight of the liquid displaced from 
the cavity to the momentum of the jet. It was at- 
tempted to derive scaling laws from these investiga- 
tions to predict penetration depth of full-scale rocket 
motors of space vehicles launched from platforms built 
over a water surface. Splash heights were also ana- 
lyzed and correlated. 


TECHNICAL MEMORANDUM X-53215 
March 12, 1965 

LAUNCH VEHICLE GUIDANCE EQUATIONS FOR 
THE SATURN IB SA-201 MISSION 

By 

Isaac E. Smith and F. Don Cooper 

George C. Marshall Space Flight Center 
Huntsville, Alabama 


ABSTRACT 


This report defines and describes the Launch Ve- 
hicle Guidance Equations for the Saturn IB SA-201 
mission, to satisfy document #15 of the Apollo Mission 
Planning and Trajectory Documentation Plan. 


TECHNICAL MEMORANDUM X-53216 
March 12, 1965 

THE SOLAR FLARE ENVIRONMENT 

By 

W. T. Roberts 

George C, Marshall Space Flight Center 
Huntsville, Alabama 


ABSTRACT 


Of all the natural environments of extraterrestrial 
space, probably the most hazardous is the radiation 


to be encountered. The three main sources of radia- 
tion are galactic cosmic rays, Van Allen belt par- 
ticles, and solar flare particles. Galactic cosmic 
rays are very high energy particles, but since their 
flux is low, they will not be considered in this analysis. 
The Van Allen belt particles will not be considered 
since we are specifically interested in the environment 
encountered at some distance from earth greater than 
about fifteen earth radii. Thus, we are left with the 
solar flare as the main contributing source of intensive 
radiation. 

Discussed is the solar atmosphere including the 
photosphere, chromosphere, and corona. Most of the 
main structures found within these areas are men- 
tioned, and some are elaborated on. 

Next are described the various iypes of solar 
flare s which have been observed . The problem of con- 
structing a solar flare "model" is discussed, and per- 
centile curves based on presently available data are 
suggested as a solution. These curves are drawn for 
the 99.9 percentile, 90 percentile, 67 percentile, and 
50 percentile groups. 

Probabilities of solar flare occurrence, and 
methods now being used as a basis for studies on solar 
flare prediction are considered. It is concluded that 
even accurate solar flare predictions cannot alleviate 
the amount of shielding required for a long duration 
extraterrestrial mission. 

A few of the complications involved in solar flare 
shielding calculations are covered, and a possible 
means of approach to these problems is suggested for 
future extraterrestrial missions. 


TECHNICAL MEMORANDUM X-53222 
March 22, 1965 

MONTE CARLO PERFORMANCE ANALYSIS COM- 
PUTER PROGRAM WITH PROGRAMMED 
MIXTURE RATIO 

By 

Roy C. Lester 


ABSTRACT 


This report presents a detailed description of a 
Fortran IV computer program which uses the Monte 
Carlo technique to establish the stage residual and a 



total payload distribution for a vehicle with up to ten 
stages. THe Monte Carlo technique employs the sam- 
pling of the independent engine and vehicle parameters, 
and, by the use of a large number of cases, is used to 
establish the probability distribution for stage residual 
and total payload. 


TECHNICAL MEMORANDUM X-53225 
March 25, 1965 

THE FIRST STAGE TRAJECTORY OPTIMIZATION 
OF THE SATURN V VEHICLE 

By 

Gary McDaniel 

George C. Marshall Space Flight Center 
Huntsville, Alabama 


ABSTRACT 


The effects of introducing optimum angle of attack 
programs into a family of Saturn V first stage trajec- 
tories are presented. These angle of attack programs 
were determined by a calculus of variations technique 
[i] such that the cutoff weight into a 100 nautical mile 
circular orbit was maximized. In tk‘.s family of first 
stage trajectories, the vehicle was constrained to fly 
a nonlifting trajectory profile (except during the time 
when small angle of attack programs were induced 
shortly after launch to initiate tilting) until the optimi- 
zation was introduced at 48, 72, and 88 seconds. The 
results of the study indicate that between 1500 and 5000 
pounds of cutoff weight in a 100 nautical mile circular 
parking orbit can be gained as compared to the optimi- 
zed nonlifting trajectory. The resulting trajectories 
show significant increases in the dynamic pressure and 
large ql«l products which led to using a q|ot| con- 
straint to reduce the structural loadings. Additional 
studies are required to determine the structural load- 
ings and the control requirements associated with the 
introduction of the optimization. 


TECHNICAL MEMORANDUM X-53235 
April 9, 1965 

SA-9 FLIGHT TEST DATA REPORT 


George C. Marshall Space Flight Center 
Huntsville, Alabama 


ABSTRACT 


This report is a presentation of certain flight me - 
chanical data obtained from the SA-9 flight test. Digi- 
tized data are presented in graphical form. Also 
included are schematic drawings showing the instru- 
ment location on the vehicle. 

The intention of this report is to present the digi- 
tized data in an easy -to-read form for use by design 
and technical personnel. This report is to supplement 
the Saturn SA-9 Flight Evaluation Report and many 
other reports published by the various laboratories. 


TECHNICAL MEMORANDUM X-52242 
April 15, 1965 

SA-201 LAUNCH VEHICLE REFERENCE 
TRAJECTORY (U) 

By 

Joseph W. Cremin and William M. Gillis 

George C, Marshall Space Flight Center 
Huntsville, Alabama 


ABSTRACT 


Presented is the official launch vehicle reference 
trajectory for the first R&D flight of the Saturn IB: 
The various launch vehicle (S-IB and S-IVB stage), 
spacecraft, and system constraints that have produced 
this trajectory shaping are reported in addition to an 
unofficial representation of the Apollo Command and 
Service Module (CSM) powered flight. The nominal 
launch vehicle mission is the injection of a partially 
fueled CSM into an earth intersecting ellipse. The SM 
is ignited after a specified codst phase and propels 
the CM to atmospheric reentry conditions, that will 
provide an evaluation of the heat shield ablator per- 
formance at high heat rate conditions. 


TECHNICAL MEMORANDUM X-53244 


H. J. Weichel 


April 15, 1965 



HARDWARE AND GUIDANCE SCHEME ERROR 
ANALYSIS FOR SA-201 EVALUATED AT S-IVB 
CUTOFF 

By 

A. W. Deaton 


ABSTRACT 


This technical memorandum presents the inac- 
curacies in velocity and displacement of the Saturn IB 
SA-201 vehicle at S-IVB stage cutoff caused by hard- 
ware and guidance scheme errors. The hardware 
errors include such effects as stabilized platform drift 
and accelerometer measurement inaccuracies. 


TECHNICAL MEMORANDUM X-53251 
April 30, 1965 

SATURN SA -9/PEGASUS A POSTFLIGHT 
TRAJECTORY 

By 

Jonathan B. Haussler and Robert H. Benson 

George C. Marshall Space Flight Center 
Huntsville, Alabama 


ABSTRACT 


This report presents the postflight trajectory for 
the Saturn SA-9/Pegasus A test flight. Fourth of the 
Block II series, SA-9 was the first vehicle to carry a 
Pegasus payload. Trajectory dependent parameters 
are given in earth-fixed, space-fixed ephemeris, and 
geographic coordinate systems. A complete time his- 
tory of the powered flight trajectory is presented at 
1. 0 sec intervals from first motion to S-I/S-IV sepa- 
ration and at 5.0 sec intervals from S-I/S-IV separa- 
tion to insertion. Tables of insertion conditions and 
various orbital parameters are included in a discussion 
of the orbital portion of flight. 


TECHNICAL MEMORANDUM X-53255 
May 5, 1965 


MEASUREMENT RELATIONSHIPS FOR PEGASUS- 
TYPE METEOROID EXPERIMENTS 

By 

Charles C. Dalton 

George C. Marshall Space Flight Center 
Huntsville, Alabama 


ABSTRACT 


In Pegasus-type meteoroid measurement satel- 
lite experiments, where vehicle attitude information 
is available without information about the meteoroid 
closing velocity vector, the probable radiant and the 
puncture flux enhancement factor for sporadic meteors 
should be based on the assumption that velocity vectors 
for ambient meteoroids are isotropically distributed. 
Then, although half of the meteoroids which impact 
onto a flat transducer are more than 45 degrees from 
the normal, those within 45 degrees puncture three 
times more effectively than the others, and half of the 
puncturing meteoroids are less than 33.5 degrees 
from the normal when the earth is below the surface 
horizon. This distribution also implies that the flux 
of puncturing meteoroids is five times greater than 
the flux of incident meteoroids with mass equal to or 
greater than the mass of the almost-puncturing mete- 
oroid with the average values of density, velocity, 
and impact angle. 

The formulas for the celestial coordinates of the 
probable radiant of a puncturing meteoroid have been 
corrected for earth-shielding. This correction is im- 
portant near the earth even under nominal oiroum- 
stances; e.g. .when the orbit height is 500 kilometers, 
the radius of the earth's disk (as seen from the ve- 
hicle) is 70 degrees. In that case, when the elevations 
of the earth's center are 0, 15, 30, 45, and 60 de- 
grees, the probable radiants are defleoted 5, 10, 20, 
39, and 80 degrees, respectively, by the earth- 
shielding effect. 


TECHNICAL MEMORANDUM X-53257 
May 11, 1965 

HEAT TRANSFER IN THE VICINITY OF A TWO- 
DIMENSIONAL PROTUBERANCE 

By 
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Ed Murply 


TECHNICAL MEMORANDUM X-53259 


George C, Marshall Space Flight Center 
Huntsville, Alabama 


ABSTRACT 


Empirical equations are presented which will es- 
timate the ratio of protuberance-to-flat-plate heat 
transfer coefficients for supersonic flow in the vicinity 
of two-dimensional protuberances in a turbulent, 
supersonic boundary layer. The equations were de- 
rived from experimental heat transfer data taken 
around several two-dimensional protuberances, i.e. , 
a 1 x 2-inch (2.5 x 5 cm) stringer, a 2 x 4-inch (5 x 
10 cm) stringer with and without a 30-degree wedge, 
and a one -quarter cylinder forebody. Important pa- 
rameters in the equations are boundary layer thickness , 
protuberance geometry, Reynolds number, and Mach 
number. Using the equations, the average deviation 
between estimated and experimental ratios is from -9 
to +15 percent. 


TECHNICAL MEMORANDUM X-53258 
May 12, 1965 

A THREE-DIMENSIONAL ASCENDING ITERATIVE 
GUIDANCE MODE 

By 

Isaac E. Smith 

George C. Marshall Space Flight Center 
Huntsville, Alabama 


ABSTRACT 


This report discusses a novel computational pro- 
cedure for implementing three-dimensional, path 
adaptive guidanoe in upper stage vaouum flight. The 
procedure is characterized by a minimum number of 
preset values, on-board computation of guidance pa- 
rameters, and reasonably close adherence to trajec- 
tories determined by calculus of variations methods, 
both as to trajeotory shape and payload capability. 
Using average gravity magnitude and direction approx- 
imations, a partially dosed form solution to the equa- 
tions of motion is obtained. The solution yields 
a "time-to-go" and a steering funotion. The steering 
funotion oonsists of a thrust attitude and a turning rate. 


May 12, 1965 

NASA’S 150-METER METEOROLOGICAL TOWER 
LOCATED AT CAPE KENNEDY, FLORIDA 

By 

John W. Kaufman 

George C. Marshall Space Flight Center 
Huntsville, Alabama 

and 

Lester F, Keene 
John F, Kennedy Space Center 


ABSTP4.CT 


This report describes the basic functional oper- 
ations of the 150-meter NASA Meteorological Tower 
Facility located in the Merritt Island Launch Area, 
Cape Kennedy, Florida. Meteorological instrumenta- 
tion and data acquisition schemes for obtaining lower 
atmospheric profile measurements of temperature, 
winds, and humidity are discussed. Specific infor- 
mation on such matters as meteorological sor.jor re- 
sponse, electronic circuitry of the meteorological 
facility, and explicit details of recorder systems has 
purposely been omitted. 


TECHNICAL MEMORANDUM X-53261 
May 14, 1965 

SOME PRACTICAL ACCURACY CONSIDERATIONS 
OF SMOKE TRAIL WIND PROFILE DATA 


By 

Dennis W. Camp 
and 

James R. Scoggins 

George C. Marshall Space Flight Center 
Huntsville, .Alabama 





ABSTRACT 


Wind velocity profiles measured simultaneously 
by the smoke trail/photographic method by different 
camera pairs, and photographs reduced twice for the 
same trail are aralyzed to determine the relative ac- 
curacy of the smoke trail wind data. Arithmetic 
means, standard deviations, and extremes are given 
for the differences between five pairs of profiles. Also 
presented are relative RMS errors in wind speeds, 
relative RMS errors in shears over height increments 
of 25, 50, 100 and 200 m, and maximum shear dif- 
ferences for each case. The computed RMS differ- 
ences in wind speeds are generally less than 1 m/sec, 
errors in wind shear vary to approximately 60 percent 
pf the values used in vehicle design studies, and rela- 
tive RMS errors in wind speeds range from approxi- 
mately 0. 1 to 0,7 m/sec. The sources of the errors 
are not isolated. These errors will probably be re- 
duced as more experience is gained in data reduction 
and quality control procedures are improved. 


TECHNICAL MEMORANDUM X-53262 
May 17, 1965 

SA-8 OPERATIONAL TRAJECTORY 

By 

Pamelia B. Pack 
and 

Gordon W, Solmon 

George C. Marshall Space Flight Center 
Huntsville, Alabama 


ABSTRACT 


This report presents the operational predicted 
trajectory for Saturn I vehicle SA-8 to be flown over 
the Eastern Test Range. Included is a discussion 
of the vehicle and mission objectives, the trajectory 
shaping and coRotraints and a brief description of the 
vehicle. Wind limit criteria and range safety aspects 
are also summarized, 

A. successful flight will insert the depleted S-IV 
stage and payload into an orbit with perigee and apo- 
gee altitudes of 510 km and 754 km, respectively. 


This orbit has a nominal lifetime of 1200 days. 

The payload consists of an Apollo boilerplate 
(BP-26) and a Meteoroid Technology Satellite (Peg- 
asus) . BP -26 is used to simulate the characteristics 
of an Apollo spacecraft whose ultimate mission will be 
a manned lunar landing and return to earth. After in- 
sertion into orbit, Pegasus will be exposed and, after 
wing deployment, will transmit micrometeoroid data 
to earth upon telemetered command. 


TECHNICAL MEMORANDUM X-53263 
May 17, 1965 

SA-9, 8 AND 10 DISPERSION ANALYSIS 

By 

Gerald Wittenstein and Jerry D, Weiler 

George C. Marshall Space Flight Center 
Huntsville, Alabama 


ABSTRACT 


This report replaces the Part II portion ' the 
SA-9 Final Trajectory. The dataused to genera* ± the 
analysis presented here are based on final SA-9 data. 
The only differences in the flights of SA-9, 8, and 10 
will be a shift in the nominal similar to the shift de- 
scribed in this report. 

A detailed discussion of dispersions, their phi- 
losophy, and application with respect to a nominal tra- 
jectory is presented, A brief description of the SA-9 
nominal trajectory is contained as a reference in this 
report. 

It is felt by the authors that the group of trajec- 
tories presented in the dispersion analysis is repre- 
sentative and yields a satisfactory envelope for the 
remaining Block H flight profiles (SA-8 and SA-10). 
Considering the vehicle ± 2 cperturbations, all mis- 
sion objectives and requirements will be met with a 
high level of confidence. 


TECHNICAL MEMORANDUM X-53264 
May 21, 1965 

LAUNCH VEHICLE CONTROL SYSTEMS RESEARCH 
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By 

J. A. Lovingood 

George C. Marshall Space Flight Center 
Huntsville, Alabama 

ABSTRACT 


be judged from individual values, but only from the 
distribution of the values which it will assume in the 
long run, i.e. , from its sampling distribution. Some 
estimate of the variance of these sample estimates is 
needed. Values of auxiliary functions required to ob- 
tain, from doubly truncated normal samples, maxi- 
mum likelihood estimates of parameters of the parent 
population and the asymptotic (large sample) vari- 
ances and covariance of these estimates are given. 


Three areas of launch vehicle control research 
are discussed: adaptive bending stabilization, adaptive 
load relief, and a minimax control problem. For each 
of these problems, factors which motivate the re- 
search, specific systems which have been investigated 
by contractors, and the status of these investigations 
are presented. In conclusion the future direction of 
the Aero-Astrodynamics Laboratory's research pro- 
gram is outlined. 


TECHNICAL MEMORANDUM X-53221 


TECHNICAL MEMORANDUM X-53266 
May 24, 1965 

LAUNCH VEHICLE GUIDANCE EQUATIONS FOR 
THE SATURN IB SA-202 MISSION 

By 

Isaac E. Smith and F. Don Cooper 

George C. Marshall Space Flight Center 
Huntsville, Alabama 


May 21, 1965 

AUXILIARY ESTIMATING FUNCTIONS FOR ABSTRACT 

DOUBLY TRUNCATED NORMAL SAMPLES 


By 

J. David Lifsey 

George C. Marshall Space Flight Center 
Huntsville, Alabama 


This report defines and describes the Launch 
Vehicle Guidance Equations for the Saturn IB SA-202 
mission, to satisfy document #15 of the Apollo Mission 
Planning and Trajectory Documentation Plan. 


TECHNICAL MEMORANDUM X-53268 


ABSTRACT 


May 25, 1965 


When sampling procedures on a random variable 
X are such that the resulting sample consists of N 
measured observations for which X^ < X < X^, and no 

information is known for X < X and X > X„, the sam- 

<x p 

pie is said to be doubly truncated at the known termi- 
nals X and X. To calculate maximum likelihood 
a p 

estimates of the mean and standard deviation of a 
normally distributed population from doubly truncated 
samples, it is necessary to solve simultaneously a 
pair of rather complex nonlinear estimating equations. 
Since every estimate is a function of the sample values 
and must be regarded as an observed value of a cer- 
tain random variable, there are no means of predicting 
in a given case , the true population value assumed by 
the estimate. The "goodness" of an estimate cannot 


A METHOD OF IMPLEMENTING CUTOFF 
CONDITIONS FOR SATURN V LUNAR MISSIONS 
OUT OF EARTH PARKING ORBIT ASSUMING A 
CONTINUOUS GROUND LAUNCH WINDOW 

By 

F. Don Cooper 

George C. Marshall Space Flight Center 
Huntsville, Alabama 


ABSTRACT 


This report demonstrates a method of imple- 
menting cutoff conditions for the Saturn V lunar 


i 17 


mission out of earth parking orbit assuming a ground 
launch window of four and one -half hours. 


TECHNICAL MEMORANDUM X-53245 
May 26, 1965 

ESTIMATION IN MIXTURES OF POISSON AND 
MIXTURES OF EXPONENTIAL DISTRIBUTIONS 

By 

A, Clifford Cohen, Jr.* 

George C. Marshall Space Flight Center 
Huntsville, Alabama 


ABSTRACT 


In the analysis of experimental data, many of the 
distributions encountered are the result of combining 
two or more separate component distributions. Es- 
timation in these compound or mixed distributions is 
therefore of particular interest to aerospace scien- 
tists. Estimators are derived for the parameters of 
a compound Poisson distribution with probability den- 
sity function 

-p x -X x 

f(x) = a — — ^ — + (i - a) , x=0, 1 , 2 ,. . . 

X.' X'. 

and for a compound exponential distribution with prob- 
ability density function 

f(x) = a(i/n)e x ^ tx + (i - a) (1/X)e X/ ^, xs 0 

where a is the proportionality factor (0 ^ a < l) and 
where ju and X are component parameters. In addition 
to the more general case in which all parameters must 
be estimated from sample data, several special cases 
are considered in which one or more of the parameters 
are known in advance of sampling. 


* Professor of Mathematics, University of Georgia, 
Athens, Georgia. Theresearchreported in this paper 
was performed under NASA Contract NAS8-11175 with 
the Aerospace Environment Office, Aero-Astrody- 


TECHNICAL MEMORANDUM X-532,71 
May 26, 1965 

ANALYSIS OF WIND TUNNEL DATA FOR SEVERAL 
BECKMAN & WHITLEY SERIES 50 AND CLIMET 
MODEL Cl-14 ANEMOMETERS 

By 

Dennis W, Camp 

George C. Marshall Space Flight Center 
Huntsville, Alabama 

ABSTRACT 


This report presents results of wind tunnel test 
data for several Beckman & Whitley Series 50 and 
Climet Model Ci-14anemometers. The primary pur- 
pose of the wind tunnel testing was the determination 
of the followingtwo response parameters: the distance 
constant and the damping ratio of the anemometers. 

The dependence of the distance constant on the 
electronic filter used to filter the output signal of the 
wind speed sensor is discussed. The importance of 
the distance constant in measuring the fluctuations of 
the wind is illustrated by showing how this parameter 
affects the frequency response curves of the anemo- 
meters. The average values of the distance constants 
were found to be 1. 14 m for the Beckman & Whitley 
Series 50 and 0.73 m for the Climet Model Ci-14, 

The damping ratio, in agreement with theory, is 
shown to be independent of wind speed. The values 
obtained for this parameter were 0. 59 for the Beckman 
& Whitley Series 50 and 0.47 for the Climet Model 
Ci-14 anemometers. 


namics Laboratory, Marshall Space Flight Center, 
Huntsville, Alabama. Mr. O. E. Smith and Mr, J. D. 
Lifsey are the NASA contract monitors. 
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TECHNICAL MEMORANDUM X-53272 
May 26, 1965 

THE EFFECT OF TRUNCATION ON TESTS OF 
HYPOTHESES FOR NORMAL POPULATIONS'' 

By 

Britain J. Williams 

George C. Marshall Space Flight Center 
Huntsville, Alabama 


ABSTRACT 


The loss of power involved in using the usual tests 
of hypotheses for complete normal populations when, 
in fact, the population is a truncated normal one is 
examined. The technique employed is to obtain an 
asymptotic expansion of the distribution of sums of 
samples of size n drawn from a truncated normal pop- 
ulation. An electronic computer is then employed to 
consider enough terms of this expansion to obtain the 
desired accuracy for the loss function, which is tabu- 
lated as a function of selected truncation points and 
sample sizes. 


* The research reported in this paper was submitted 
as a Ph. D dissertation in statistics at the University 
of Georgia, Athens, Georgia. This research was 
performed under NASA Contract NAS8-11175 with the 
Aerospace Environment Office, Aero-Astrodynamics 
Laboratory, Marshall Space Flight Center, Huntsville, 
Alabama, Mr, 0, E, Smith and Mr, J. D, Lifseyare 
the NASA contract monitors. The author is currently 
with the RCA Service Company, Missile TestProiect, 
Patrick AFB, Florida. 


TECHNICAL MEMORANDUM X-63273 
May 27, 1965 

SPACE ENVIRONMENT CRITERIA GUIDELINES 
FOR USE IN SPACE VEHICLE DEVELOPMENT 
( 1965 Revision) 

By 

Robert E, Smith 


George C. Marshall Space Flight Center 
Huntsville, Alabama 

ABSTRACT * 


This document provides guidelines on interplane- 
tary space, terrestrial space, near -Venus space, 
near -Mars space, lunar atmosphere and surface, 
Venus atmosphere and surface, and Mars atmosphere 
and surface environmental data applicable for Marshall 
Space Flight Center space vehicle development pro- 
grams and studies related to future NASA programs. 

This report establishes design guideline values 
for the following environment parameters: (1) me- 
teoroids; (2) secondary ejecta; (3) radiation; (4) 
gas properties; (5) magnetic fields; (6) solar radio 
noise; (7) winds; (8) wind shear; (9) clouds; 
(10) ionosphere; (11) albedo; (12) planetary sur- 
face conditions; (13) planetary satellites; (14) com- 
position of the planetary atmospheres; (15) temper- 
ature; and (16) astrodynamic constants. Additional 
information may be located in the references cited. 

Extensive use was made of the data prepared by 
the Planetary Atmosphere Section, ASTD, Space En- 
vironment, Manned Spacecraft Center, Houston, Texas, 
to insure compatibility of both development and study 
effort between Marshall Space Flight Center and the 
Manned Spacecraft Center, especially in those areas 
where there are insufficient data to make definite 
conclusions. 


* This document supersedes and should be used in lieu 
of NASA TMX-53142, September 30, 1964. 


TECHNICAL MEMORANDUM X-53275 
June 8, 1965 

PRELIMINARY SURVEY OF OPTIMUM THRUST 
LEVELS FOR THE HG-3 HIGH PRESSURE 
LOX/HYDROGEN ENGINE AS APPLIED TO 
AN UPRATED SATURN V VEHICLE 

§7 

Emsley T. Deaton, Jr. & Neva H. Huffaker 

George C. Marshall Space Flight Center 
Huntsville, Alabama 
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ABSTRACT 


Results of a design point survey are presented, 
determining the optimum thrust levels for the HG-3 
engine when used for mission profiles applicable to the 
Saturn V class of vehicles. The configuration desig- 
nated MLV-SAT V-3 from the Saturn V Improvement 
Studies was used for the trajectory analysis, with the 
J-2 engines replaced by the HG-3 engines in the second 
and third stages. 

A new design point engine was selected with the 
results of the survey , and trajectory optimizations for 
appropriate mission profiles were completed for this 
new engine. 


TECHNICAL MEMORANDUM X-53283 
June 18, 1965 

STATIC LONGITUDINAL STABILITY AND AXIAL 
FORCE CHARACTERISTICS OF A 4. 33-CALIBER 
BLUNT CYLINDER WITH A 30° SLOTTED 
FRUSTUM AFTERBODY OVER A MACH 
NUMBER RANGE OF 0. 50 TO 4. 96 

By 

Robert E. Pitcock 

George C. Marshall Space Flight Center 
Huntsville, Alabama 


ABSTRACT 


Wind tunnel tests were conducted to determine the 
static longitudinal stability and axial force character- 
istics of a proposed Saturn IB camera capsule config- 
uration consisting of a blunt cylinder with a 30° slotted 
frustum afterbody. Tests indicated this configuration 
to be statically stable about a center of gravity located 
2,43 calibers forward of the base throughout the test 
Mach number range of 0.50 to 4.96. Maximum sta- 
bility margin exists at approximately Mach 0,95 with 
minimum stability margin existing from about Mach 
number 1,20 to 4,96, 


TECHNICAL MEMORANDUM X-53288 
June. 28, 1965 


DETERMINATION OF AERODYNAMIC FOJtCE AND 
HEAT TRANSFER PROPERTIES FOR A CONCAVE? 
HEMISPHERICAL SURFACE IN FREE 
MOLECULE FLOW 

By 

Clarence Ray Wimberly 

George C. Marshall Space Flight Center 
Huntsville, Alabama 


ABSTRACT 


The methods and equations are presented for the 
aerodynamic force and heat transfer properties of a 
concave hemispherical surface at an arbitrary angle 
of attack in free molecule flow. The effect ofmulti- 
ple reflections is taken into account by assuming a 
cosine distribution after each collision. It is assumed 
that the emission of molecules from the surface is 
purely diffuse , and that the reflecting molecules are 
perfectly accommodated to the surface conditions. The 
resulting equations are not in closed form, but have 
been solved numerically for some typical values of the 
molecular speed ratio and surface-to-ambient-tem- 
perature ratio. It is found that under the above con- 
ditions the heat transfer characteristics are the same 
as those of a convex hemisphere. The drag is in- 
creased slightly over its counterpart, while the total 
lift is shown to be very small and would be zero except 
for interreflected molecules. 


TECHNICAL MEMORANDUM X-53289 
June 29, 1965 

HARDWARE AND GUIDANCE SCHEME ERROR 
ANALYSIS FOR SA-202 EVALUATED AT 
S-IVB CUTOFF 

By 

P. D. Brandon 


ABSTRACT 


This technical memorandum presents the inaccu- 
racies In velocity and displacement of the Saturn IB 
SA-202 vehicle at S-IVB stage cutoff oaused by hard- 
ware and guidance scheme errors. The hardware er- 
rors include such effects as stabilized platform drift 
ana accelerometer measurement inaccuracies. 
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B. PRESENTATIONS 

SPHERICAL BALLOON WIND SENSOR BEHAVIOR 

By 

James R. Scoggins 


ABSTRACT 


Response characteristics of freely rising super- 
pressured spheres of different configurations are dis- 
cussed. Wind profile data measured by the super- 
pressured balloon method, the smoke trail method, the 
AN/GMD-1 rawinsonde system, and from low-level 
oper-air tests were used in the analysis. Results are 
reported on data measured at Huntsville, Alabama, at 
night during stable conditions to an altitude of 120 m 
and at Cape Kennedy to an altitude of 12 km. The re- 
sults show that: (1) the average drag curve for a 
freely rising 2-m diameter, smooth superpressured 
sphere differs considerably from the drag curve ob- 
tained in wind tunnels using smaller spheres; the av- 
erage value of the drag coefficient is larger over all 
Reynolds numbers except near the transition region; 

(2) the drag coefficient for roughened spheres is 
nearly independent of the Reynolds number but de- 
creases slightly as the Reynolds number decreases; 

(3) the addition of surface roughness elements re- 
duces the aerodynamically induced horizontal motions 
of the smooth sphere; and (4) the average aerody- 
namic lift force, which acts primarily in the horizontal 
direction and is responsible for the aerodynamically 
induced horizontal motions, is negligible indicating 
that it does not act in any preferred direction. 


SUMMARY OF SATURN I BASE THERMAL 
ENVIRONMENT 

By 

Robert G. Payne and Ira P. Jones, Jr. 


ABSTRACT 


Research and development flights of the Saturn I 
have established the base thermal environment of 
eight-engine LOX/RP-1 propelled vehicles. Plume 
radiation levels, total heating rates, and base gas 
temperatures have been measured on the first stage of 
the Block I and II vehicles. Also some heating data 
were obtained for the six-engine LOX/LH z , S-IV 


second-stage base. Discussion of flight instrumenta- 
tion, data reduction techniques , and evaluation prob- 
lems is presented along with a correlation of the flight 
results to both long and short duration model tests. 
Maximum base heating rates occurred at liftoff for 
both total and radiation heatingwith a lower secondary 
peak occur ring when the reverse flow from the inboard 
engines became choked. Maximum flame shield total 
heating rates were approximately three times (75-100 
W/cm 2 ) the heat shield values for the S-I base. Base 
gas temperatures varied widely with maximum values 
of 1300°K for the heat shield and 2100 °K for the flame 
shield. Block n flight data for the S-I stage were gen- 
erally lower than the Block I results. S-IV stage base 
heat rates varied with radial distance from the center 
and flight time from 2 to 0.5 W/cm 2 and agreed well 
with short duration model tests results. 

Orbital Tracking and Decay Analysis of the Saturn 
I Flights. H. F. Kurtz, Jr. , Chief, Operations 
Studies Branch, M. Naumcheff, Chief, Orbit Deter- 
mination Section, and A. R. McNair, Chief, Mission 
Studies Section, Marshall Space Flight Center, 
Huntsville, Alabama. 

The successful skin -tracking of the first three 
Saturn I orbital flights by a network of FPS-16 and 
FPQ-6 radars is discussed, and the effective atmos- 
pheric drag force acting on the orbiting vehicles is 
presented as determined through use of the radar data 
for orbit determination. The feasibility of skin- 
tracking the Saturn vehicle in orbit was proven on the 
flightof SA-5, the first Saturn I orbital vehicle. Skin- 
track data has since reliably been obtained on vehicles 
SA-6 and SA-7. Error characteristics of the data are 
presented, which show it to be of comparable quality 
to that obtained with beacon tracking by the same ra- 
dars. The. root -mean-square of the gross error resi- 
duals obtained during orbit determinations using the 
skin -track data averaged about 30 m in slant range and 
0. 03 degrees in the azimuth and elevation angles. 
Atmospheric drag data derived from the orbit deter- 
minations are presented, with a discussion of the 
methods used. The data show current atmospheric 
density values and drag characteristics of the Saturn 
vehicles, and lead to an increase in the accuracy and 
confidence of orbital lifetime predictions for such ve- 
hicles. 


Presented to AIAA, Huntsville, Alabama, Feb- 
ruary 15 - 17, 1965, 
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ON THE CALCULATION OF DYNAMIC STABILITY 
DERIVATIVES FOR POINTED BODIES OF 
REVOLUTION IN SUPERSONIC FLOW 

By 

Max F. Platzer 

NASA, George C. Marshall Space Flight Center 
and 

Gilbert H. Hoffman 
Lockheed Missiles and Space Company 
Huntsville, Alabama 


ABSTRACT 


This paper reviews two linearized approaches to 
calculate the dynamic stability derivatives of slowly 
oscillating pointed bodies of revolution in supersonic 
flow, namely, (a) a generalized Karman-Moore ap- 
proach and (b) a generalized Adams-Sears approach. 

The Karman-Moore approach is a purely numeri- 
cal method, whereas Adams-Sears' quasi-slender body 
theoiy leads to analytical closed form solutions for 
the stability derivatives. Four configurations are in- 
vestigated in more detail, the cone, the convex and 
concave parabolic ogive and a typical reentry body. 
Comparison with other theories and experiment is 
made, thus giving an insight into the range of appli- 
cability of both theories. 


Presented at Second Technical Workshop on Dy- 
namic Stability Testing, Arnold Engineering Develop- 
ment Center, Tullahoma, Tennessee, April 21 - 23, 
1965. 


ON LINEARIZED COMPRESSIBLE FLOW OVER 
PULSATING AND OSCILLATING BODIES 

By 

M, F. Platzer 


ABSTRACT 


The paper presents a new and unified approach to 
the calculation of linearized sub- and supersonic flow 
over pulsating and oscillating low aspect ratio bodies 
extending K. Oswatitsch's and F, Keune's theories 


for steady flow to these unsteady flow cases. 

It is shown that the sub- and supersonic flow over 
low aspect ratio wings at zero angle of attack, whose 
skin executes symmetric pulsations, is in a first ap- 
proximation given by two terms; namely, a crossflow 
and a spatial influence. This spatial influence reduces 
the flow around pulsating wings to the flow around the 
equivalent pulsating body of revolution. A similarly 
characteristic structure of the flow field is found also 
for the higher order flow terms for both pulsating and 
oscillating flow. To examine its range of validity this 
approximation theory is applied to configurations for 
which exact solutions of the unsteady linearized po- 
tential equation can be found. Such exact solutions do 
exist for the infinitely long tube and ribbon exhibiting 
a harmonic standing pulsation or oscillation and the 
paper discusses these solutions in more detail. Fur- 
thermore, the slowly oscillating body of revolution is 
studied by investigating three different body configu- 
rations; namely, a cone, a convex parabolic ogive and 
a concave parabolic ogive. The stability derivatives 
are evaluated and compared with results from other 
theories. 


Presented at the 7th International Aeronautical 
Congress, Paris, June 14 - 16, 1965. 


SOME RESULTS OF AN ANALYTICAL AND 
EXPERIMENTAL INVESTIGATION TO 
DETERMINE HEAT TRANSFER IN THE 
NOZZLE OF SOLID PROPELLANT 
ROCKET MOTORS 

By 

Ewell Roger Heatherly 


ABSTRACT 


The combustion chamber of a solid propellant 
rocket motor contains a mixture of fuel and oxidizer 
in solid form arranged in a pattern to produce a nearly 
constant chamber pressure and temperature during 
the firing period. The firing duration is a function of 
amounts and types of fuel, oxidizer, and additives, as 
well as grain size and configuration of the combustion 
chamber and the nozzle. 

The nozzle controls the quantity of gas flow and 
guides the thrust. Since the gas temperature is above 
the values for which metals retain their structural in- 
tegrity, some provision must be made to cool the 
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nozzle. Irt general, there are three basic methods of 
cooling nozzles which are functions of the firing du- 
ration. For a very short duration of firing, the heat 
sink method of thermal protection is used , talcing ad- 
vantage of transient heat transfer. For a moderate 
length of firing, an ablative nozzle may be used. The 
heat is carried away by actual sublimation, melting, 
evaporating, or charring of the surface of the nozzle. 
The ablation products tend to act as a buffer layer for 
areas downstream, similar to transpiration cooling. 
For long firing durations , transpiration cooling may 
be used. 

The design of solid propellant motors is dependent 
on an accurate prediction of the nozzle wall heat flux. 
The heat flux into the nozzle wall or surface of the pro- 
tective layer is a function in various degrees of the 
following: 

1. Physical properties of gas constituents. 

2. Chemical kinetics during expansion. 

3. Metal additives in the propellant. 

4. Local flow pattern and properties in the 
nozzle. 

5. Boundary layer thickness and growth. 

The prediction of the heat transfer to the heat sink 
nozzle wall is especially complicated in that it re- 
quires some knowledge of the accumulation of metal 
oxide on the walls (aluminum, in this investigation), 
which increases resistance to heat flow to the non- 
ablative nozzle. 

hi the present study, which employs both experi- 
mental and analytical methods, a mathematical model 
is developed for the heat sink nozzle which combines 
the heat transfer due to convection, conduction, and 
radiation and also includes the heat flow resistance due 
to the oxide layer. Ten subscale solid propellant mo- 
tors were test fired to obtain experimentally the noz- 


zle temperatures which verified the theoretical tem- 
peratures computed by use of the mathematical model. 
Computer programs were used to predict nozzle tem- 
perature, heat transfer coefficient, resistance to heat 
flow due to oxide deposit, and gas properties. The 
selection of the proper technique for each mode of heat 
transfer required a literature search and an evaluation 
of the relative merits of various methods of study. 

Although at first the closed form heat conduction 
solutions already available [1,2,3] were considered, 
these were subsequently discarded because the tesLs 
showed that the boundary conditions were dependent on 
time. An available numerical computer program cap- 
able of treating a wide variety of boundary conditions 
has been outlined in Chapter II-A and was chosen for 
use in this analysis. 

Several empirical methods [4-9] were available 
for estimating the coefficient of convective heat trans- 
fer, but these were found to be insufficient to describe 
the details of the environment to which the nozzle was 
subjected. Bartz and his co-workers [10-14] solved 
the integral boundary layer momentum and energy 
equations to produce nozzle convective heat transfer 
relations. Because of its wide use in industry and 
good correlation with test data, this method was used 
in the calculations. 

The contribution of radiative heat transfer to the 
nozzle heating problem is generally small in compari- 
son to that of convective heating. References 15 
through 18 are used in the present study to calculate 
the radiative heat transfer due to the aluminum oxide 
particles and the products of combustion. 


Presented in partial fulfillment for M. S. at the 
University of Alabama, Tuscaloosa, Alabama, June 
1965. 
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